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

    End-to-end computational design for an EUV solar corona multispectral imager with stray light suppression

    2024-03-04 03:47:20JinmingGaoYueSunYinxuBianJilongPengQianYuCuifangKuangXiangzhaoWangXuLiuXiangqunCui
    天文研究與技術(shù) 2024年1期

    Jinming Gao, Yue Sun, Yinxu Bian, Jilong Peng, Qian Yu, Cuifang Kuang*,Xiangzhao Wang, Xu Liu, Xiangqun Cui*

    1State Key Laboratory of Extreme Photonics and Instrumentation, College of Optical Science and Engineering,Zhejiang University, Hangzhou 310027, China

    2ZJU-Hangzhou Global Scientific and Technological Innovation Center, Zhejiang University, Hangzhou 310027, China

    3Nanjing Institute of Astronomical Optics & Technology, Chinese Academy of Sciences, Nanjing 210042, China

    4Beijing Environmental Satellite Engineering Institute, Beijing 100094, China

    Abstract: An extreme ultraviolet solar corona multispectral imager can allow direct observation of high temperature coronal plasma, which is related to solar flares, coronal mass ejections and other significant coronal activities.This manuscript proposes a novel end-to-end computational design method for an extreme ultraviolet (EUV) solar corona multispectral imager operating at wavelengths near 100 nm, including a stray light suppression design and computational image recovery.To suppress the strong stray light from the solar disk, an outer opto-mechanical structure is designed to protect the imaging component of the system.Considering the low reflectivity (less than 70%)and strong-scattering (roughness) of existing extreme ultraviolet optical elements, the imaging component comprises only a primary mirror and a curved grating.A Lyot aperture is used to further suppress any residual stray light.Finally, a deep learning computational imaging method is used to correct the individual multi-wavelength images from the original recorded multi-slit data.In results and data, this can achieve a far-field angular resolution below 7", and spectral resolution below 0.05 nm.The field of view is ±3 R⊙ along the multi-slit moving direction, where R⊙represents the radius of the solar disk.The ratio of the corona's stray light intensity to the solar center's irradiation intensity is less than 10-6 at the circle of 1.3 R⊙.

    Keywords: EUV solar corona imager; Curved grating; Stray light suppression; Computational multispectral imaging

    1.INTRODUCTION

    Solar coronal observation is very important for the protection of both spacecraft and planet Earth.The EUV band is significantly different to visible light and other bands of electromagnetic radiation, with a shorter wavelength and higher photon energy[1-3].Because this band is strongly absorbed by Earth’s atmosphere, EUV telescope imagers work best in space.An EUV solar corona imager can enable direct observations of high temperature coronal plasma, and the most attractive for this application are multispectral imagers, consisting of a wide field imager and a spectrometer.Most materials absorb EUV radiation,preventing the use of prisms and transmissive optical elements.Additionally, due to the low EUV reflectivity of materials, it is unrealistic to add too many reflective elements[4,5]into EUV multispectral imagers.For example,a basic reflective Fourier transform imaging spectrometer,in which light is reflected several times, is inappropriate for EUV light[6-8].Some EUV thin-film filters have been used in the EUV imager on the Solar Dynamics Observatory satellite and the Atmospheric Imaging Assembly satellite, but these only have a spectral resolution of approximately 0.5 nm[9-12].

    With advancements in multilayer film and grating fabrication techniques, curved gratings have been introduced into EUV spectrometers.These curved gratings have unique diffraction capabilities, enabling spectral imaging and dispersion with only a single grating, reducing the system's reflective and diffractive surfaces significantly[13-15].Consequently, an EUV solar corona multispectral imager would combine a telescope system with a curved grating[16-18].The telescope images the incident light into a slit, which serves as both the field stop for the telescope and the slit for the spectrometer[19].The curved grating here works simultaneously as both imager and grating, removing the requirement to collimate mirrors, thus greatly improving the system sensitivity.The system can achieve a wide field-of-view (FOV) through scanning movement of the slits[20-22].The CORONA project(America’s first satellite program) employed a planar grating combined with a focusing mirror[23-25], able to perform multispectral imaging with high spectral resolution and wide FOV[26].

    For an EUV solar corona multispectral imager with high sensitivity and efficiency, we propose an end-to-end computational design method with separate optical elements and a deep learning computational image recovery algorithm.For suppressing stray light from the solar disk,the external occulter (EO) is innovatively designed, with an annular light aperture and an inclined spherical mirror with a central opening.For the imaging component, the simplest configuration is adopted, with only a primary mirror and a curved grating, which can achieve both high spatial and spectral resolution in a wide FOV.A unique multi-slit arrangement greatly improves the scanning efficiency.In addition, enhancement by deep learning (DL)methods can correct distortion caused by the curved grating.The entire system comprises only two optical reflective surfaces, ensuring high transmission efficiency.

    2.RESULTS

    2.1.The EUV Solar Corona Multispectral Imager System

    The main specifications of the EUV solar corona multispectral imager system are listed in Table 1.The EUV spectral ranges is 94 nm to 109 nm, namely 94 nm, 97.5 nm,97.7 nm, 102.5 nm, 103.2 nm, 103.7 nm and 109 nm.The half-angle of FOV along the multi-slit moving direction is 3R⊙a(bǔ)nd the half-FOV along the slit is 2R⊙.The spectral resolution is better than 0.05 nm, and the far-field angle resolution is better than 7".The design uses 5 slits for parallel-scanning spectral imaging.The design of the EUV solar corona multispectral imager is shown in Fig.1, including an EO structure, entrance pupil stop, offaxis primary mirror, multi-slit plane, spectral dispersion component, and Complementary Metal-Oxide-Semiconductor (CMOS) image sensor.As stray light beams from the solar disk are very strong, they should be prevented from entering the aperture, so an outer stray light suppression component is required before the imaging system.The primary mirror is then placed after the entrance pupil stop.The primary mirror is an off-axis parabolic mirror and is used to capture the first image of the solar corona, with the 5-slit multi-slit mask placed at the imaging plane of the primary mirror.A curved grating is used for multispectral imaging, conjugately imaging the first image plane to the CMOS image sensor.A Lyot stop is placed at the exit pupil to suppress stray light.The CMOS image sensor is placed on a tilted plane.Fig.2 shows the quantitative imaging performance of the system, including root mean square (RMS) data, spot radii, and modulation transfer functions (MTFs).In Fig.2A, B and C, the RMS spot radii are all below 10 μm, where the pixel size of the CMOS image sensor is 25 μm.Additionally, the MTF is better than 0.2 at 40 lp/mm across all FOVs, at wavelengths of 94 nm, 102.5 nm, and 109 nm.

    Fig.1.Schematic of the EUV solar corona multispectral imager.(A) Opto-mechanical structure, from the outer stray light suppression component to the CMOS image sensor.(B) Slit mask for multi-slit scanning used in EUV multispectral imaging.(C) Optical light path, including outer stray light suppression and imaging components.

    Fig.2.Image performance of the EUV solar corona multispectral imager system.(A) RMS spot radii change with wavelength under different FOVs.(B)-(C) RMS spot radii versus FOV at different wavelengths.(D)-(F) MTF of the optical system at wavelengths of 94 nm, 102.5 nm, and 109 nm.

    Table 1.EUV solar corona multispectral imager system specifications

    2.2.Stray Light Suppression

    To prevent stray light radiation from the solar disk,an outer stray light suppression shielding system is employed.The opto-mechanical structure is shown in Fig.3A, with two separated components.The front body includes a circular obscuration aperture (COA), while the rear body employs an annular spherical mirror aperture(ASMA).Fig.3B shows the irradiance distribution of the annular incident window, and Fig.3C shows the irradiance distribution on the primary mirror of the coronagraph, within the FOV of 0-2R⊙.Results indicate that rays within the 0-1.3R⊙FOV are imaged at the annular aperture and kept out of the imaging apparatus, where the ratio of the corona's stray light intensity to the solar center's irradiation intensity is less than 10-6at the circle of 1.3R⊙.

    Fig.3.The EO.(A) Mechanical structure of optics.(B) Irradiance distribution at the annular aperture.(C) Irradiance distribution at the image plane of the primary mirror and the multi-slit plane.

    There are two main sources through which stray light can be introduced into the imaging system, which can reduce the image contrast of the CMOS image sensor.One source is light diffracting and scattering at the edge of each element.Another, as the EUV wavelength is short, is from light scattering off any roughness or defects in the mirror and grating.To minimize the effects of this,a Lyot stop is placed at the exit pupil of the imaging component, as shown in Fig.4A.Calculations show that the diffracted light forms an annular illumination pattern witha diameter of 8.9 mm at the location of the Lyot stop.Consequently, the optical diameter of the Lyot stop is set to 8.7 mm to completely block the diffracted light.The Lyot stop is positioned 178.2 mm before the image plane with a tilt angle of -6.2°.The system vignetting caused by the Lyot stop is shown in Fig.4B.The maximum unvignetted light fraction of the system is over 95%, which will not influence the image quality.Finally, stray light is suppressed less than 10-6at the circle of 1.3 R⊙.

    Fig.4.Lyot stop to suppress stray light rays in the imaging system.(A) Position of the Lyot stop.(B) Vignetting calculation of the system caused by the Lyot stop.Please note that the y-axis here shows only 0.9-1.

    2.3.Primary Mirror

    A primary mirror is used for the first imaging, which is also the multi-slit plane.Here, an off-axis parabolic mirror works as the primary mirror, similarly to a telescope objective, the design parameters of which are given in Table 2.The radius of the parabolic mirror is -810 mm,the conic coefficient is -1, the decentration is 50 mm,and the clear aperture is 55.4 mm.The surface sagittal height profile and the cross-sectional sagittal height profile of the primary mirror are shown in Fig.5.

    Fig.5.Sagittal height of the primary mirror.(A) Surface sagittal height diagram.(B) Cross-sectional sagittal height diagram.

    Fig.6.Sagittal height diagram of curved grating.(A) Surface sagittal height diagram.(B) Cross-sectional sagittal height diagram.

    Table 2.Primary mirror parameters

    2.4.Curved Grating

    The parameters of the curved grating are divided into two parts: surface shape parameters and engraved line parameters.The surface shape is composed of an ellipsoid base surface superposed with XY extended polynomials.The surface shape expression is expressed as

    In Eq.(1) and Eq.(2),x,yandzare the intercept coordinates of the intersection point between the light ray and the surface, whilea,b, andcare the coefficients defining the major and minor axes of the elliptical shape.The grating line density is given by

    The specific surface parameters of the grating are shown in Table 3, Table 4, and Table 5.The surface sagittal height diagram and the cross-sectional sagittal height diagram of the grating are shown in Figure 6.

    Table 3.Grating surface parameters

    Table 4.Surface sagittal height polynomial parameters

    Table 5.Grating line density parameters

    Fig.7 displays the grating line density corresponding to the center position of the spot, indicating that the groove density varies with the relative position in the Y direction.The relative landing point position of the edge ray of the target FOV on the diffraction grating surface is-13.31 mm.The corresponding grating line density at this point is 2034 lines/mm.The relative landing point position of the lower edge ray of the target FOV on the diffraction grating surface is 13.31 mm, the grating line density corresponding to this point is 1970 lines/mm.Using Eq.(3) to integrate the spot range, the total number of slots obtained is 36 432, which meets the system's requirements for spectral resolution.

    Fig.7.Grating line density distribution map.The diameter of the curved grating is 40 mm.As examples, the point position at-13.31 mm on the diffraction grating surface has a grating line density of 2034 lines/mm, while the point position at 13.31 mm has a grating line density of 1970 lines/mm.

    2.5.Computational Multispectral Image Correction

    Each original spectral image is stitched from slit-bar data, which is directly recorded by the CMOS image sensor.These stitched images contain distortion commonly found in slit-scanning spectral imagers.A DL method with a convolutional neural network (CNN) is used to correct imaging distortion for the stitched spectral images, as shown in Fig.8.We use a total of 5 000 sets of data to train and test the performance of our computational multispectral image recovery, constructed within the PyTorch DL framework.The input and target values of the CNN are distorted corona images and ideal coronal images respectively.4 500 sets are used for training the CNN network and 500 sets are used for testing.A total of 20 epochs are performed, and the entire calculation process is performed on a desktop computer running the Microsoft Windows 10 operating system, with Core i7-7700K CPU @ 4.2GHz (Intel) with 64GB RAM and GeForce RTX2080Ti GPU (NVIDIA) hardware.The training time is approximately 2 hours, and the correction time of a single image is less than 200 ms.Test results are shown in Fig.8.In 500 sets of test data, the average structural similarity value of the corrected images and raw image reaches 0.898, which is higher than the 0.769 of the distorted images.

    Fig.8.Deep learning distortion correction for multispectral images correction.The reddish-brown area in the figure shows perfectly circular masks.(A) Imaging FOV diagram.(B) Simulated solar corona raw image.(C) Stitched image after slit scanning,from slit-bar data at 94 nm, with distortion.(D) Output of the deep learning CNN, without distortion.

    3.DISCUSSION

    There are two important considerations for an EUV solar corona multispectral imager, which can easily be overlooked.The first is to suppress stray light.Light from the solar disk is very strong, and this must be excluded from the imaging component.Additionally, to avoid internal scattering within the apparatus, the edges of all apertures should be crafted with minimal imperfections.All mechanical structures should be coated in light-absorbing material to minimize reflected light.Meanwhile, the mechanical structure should be strong and any vibrations should be sufficiently dampened to survive a rocket launch.Therefore, the mechanical structure should be designed with multiple chambers, with good inner space distribution, to give the instrument increased stability.This can also limit the propagation of scattered light within each narrow individual chamber, preventing the proliferation of stray light.

    The second consideration is to design and arrange suitable slits, as detailed in Section 5.The slit width should match the optical resolution and the pixel size of the CMOS image sensor.If the slits are too narrow, diffraction and exposure time can become problematic, but slits which are too wide will reduce spatial resolution to unacceptable levels.Additionally, the energy distribution of the EUV spectrum peaks at a specific wavelength.The arrangement of slits must be considered according the wavelength distribution produced by the curved grating.

    4.CONCLUSION

    In this paper, we propose an end-to-end computational design method for an EUV solar corona multispectral imager with high sensitivity and efficiency.It uses a design with separate optical elements and a DL computational image correction algorithm.The design uses an external EO with an annular light aperture and an inclined spherical mirror with a central opening, to suppress stray light from the solar disk, with a unique multi-slit arrangement to greatly improve scanning efficiency.The entire system comprises only two optical reflective surfaces, ensuring high transmission efficiency.Following image acquisition,distortion from the curved grating can be corrected using a CNN.From the design results, the EO enhances stray light suppression, and the spectrometer exhibits outstanding spatial and spectral imaging performance.The farfield angular resolution of this design is better than 7",and the spectral resolution is better than 0.05 nm.The FOV is ±3R⊙a(bǔ)long the multi-slit moving direction.The ratio of the corona's stray light intensity to the solar center's irradiation intensity is less than 10-6at the circle of 1.3R⊙.We believe that this novel method can lead to the design and development of further new EUV spectral imaging systems.

    5.METHODS

    5.1.End-to-end Computational Design

    An end-to-end computational design for the EUV solar corona multispectral imager is presented in Fig.9.The physical imaging model is generalized, from the object ‘ solar corona’ to the final recovered images.In Fig.9B, a stepwise approach is adopted based on imaging requirements, culminating in the final imaging criteria.To block the solar surface radiation from the coronagraph system as early as possible, the design of the EO and the vignette analysis are carried out first.Then, to obtain a clear image of the solar corona at the slit position, the telescope primary mirror is designed computationally.To reduce the full-field scanning imaging time and eliminate the overlap of imaging spectral lines in different FOV, it is necessary to determine the number of slits and slit arrangement.Finally, DL computational image recovery is used to correct aberrations and improve image quality.

    Fig.9.End-to-end computational design for the EUV solar corona multispectral imager.(A) The end-to-end roadmap from the object ‘solar corona’ to the final recovered images.(B) Each design point, including stray light suppression, telescope objective, multislit mask, curved grating and the deep learning recovery algorithm, along the end-to-end roadmap.

    5.2.Design of Stray Light Suppression Components

    A schematic diagram of the EO stray light suppression structure is shown in Fig.10.The structure consists of a COA and an ASMA with an inclination angle ψ.The center of the COA is an occluder, and the center of the ASMA is an aperture with a radius ofr0.A spherical mirror with a specific tilt angle, ψ, can image the sun at the annular optical entrance, thereby reflecting it out of the solar corona imager.The inner radius of the COA isr1, the outer radius isr2, and the positionsE0andM0are separated by distance,L.If the image of the sun is located onE0, the distancedbetween the image center and the COA center needs to satisfy the formula where,

    Fig.10.EO stray light suppression structure and design.It consists of a COA and an ASMA with an inclination angle.The center of the COA is an occluder, with an aperture at the center of the ASMA.The spherical mirror images the sun at the annular optical entrance, thereby reflecting it out of the solar corona imager.

    his the half-height of the image of the sun, so that Ψ satisfies the formula

    Here, θ0is the half-FOV angle of the solar surface image, and θ1is the maximum half-FOV angle of the solar coronal image.All light rays from 0- θ0will be intercepted by the ASMA, will be imaged on the surface of the annular entrance where the COA is located, and will be reflected out of the solar corona imager.All the solar coronal light from θ0- θ1enters the system, while the light greater than θ1will be partially intercepted by the ASMA.If the pointing accuracy of the satellite isφ, thenr1andr2satisfy the equation

    5.3.Primary Mirror Design

    To meet the desired spatial and spectral resolution requirements, clear imaging at the primary image plane of the system is needed, i.e., at the slit position.The main mirror adopts an off-axis parabolic design[27].When determining the parameters of the off-axis paraboloid, the following factors should be comprehensively considered: (1)Since a multi-slit plane is placed at the primary image plane of the primary mirror and an image of the corona region is obtained through scanning the multi-slit device,the imaging focusing performance (spot focusing diameter) of the primary mirror needs to match the width of the slit.(2) The focal length of the spectral imaging system only depends on the focal length of the primary mirror and the magnification of the curved grating.Therefore,the focal length design of the off-axis parabolic primary mirror needs to match the spatial resolution[28].(3) To avoid mechanical interference between the optical path and optical elements, there should be sufficient space for mechanical structure assembly and adjustment.The offaxis parabolic mirror needs to be sufficiently off-axis.At the same time, the off-axis paraboloid serves as the front telescope system of the spectral imager, and its design needs to meet the imaging focusing performance of the target FOV at the final image plane[29-31].

    5.4.Slit Distribution Design

    The multi-slit plane is located at the primary image plane and consists of straight slits with unequal spacing,as shown in Fig.11.From the focal length of the telescope imaging system, the field diaphragm size at the first image plane can be calculated as

    Fig.11.Schematic diagram of the multi-slit mask, where ‘px’ denotes size in pixels.The spectral dispersion distributions corresponding to the five slits are independent and are located together on the CMOS image sensor.The spectral bandwidths will partially overlap, but the spectral line positions do not overlap.

    where,YFSis the width in theYdirection andXFSis the length in theXdirection.In Fig.11,Nw+(N-1)ds,i=YFS, wherewis the slit width,ds,iis the slit spacing, andNis the number of slits.The image of the slit through the grating has a widthw1on the CMOS image sensor.The spectra produced byw1will partially overlap, limiting the resolvable spectral bandwidth.The spectral resolution corresponding to the geometric width of the slit image should satisfy

    Multiple slits enable different spectral line segments at different FOVs to be imaged on the CMOS image sensor simultaneously.To reduce the time cost, multiple slits should be involved in imaging while ensuring spectral resolution.The slit spacingds,ican be calculated with the formula

    where,ds,iis the spacing between i slits,dm,iis the minimum spacing between slitsiandi+1 on the image plane,andβis the magnification of the curved grating.

    5.5.Curved Grating Design

    To correct chromatic aberration and obtain high spectral resolution imaging, the spectral dispersion system uses a variable line spacing curved grating.The aberration correction principle diagram is shown in Fig.12.The pointA0(xa,ya,0) is the center of the slit, andAA0is any off-axis light source parallel to the grating line.The pointP(x,y,z)is any point on the grating surface.After the incident light raysAPandAOpass through the curved grating, they are imaged at pointB(xb,yb,zb).B0(xb,yb,0) is the projection of pointBon the planeXOY.The angleibetween the lightA0Oand theX- axis is the incident angle, and |OA0|=rAis the incident arm length of the grating.The angleθbetween the lightB0Oand theX-axis is the diffraction angle, and |OB0|=rBis the exit arm length of the grating.The surface sag of the grating can be expressed using Eq.(1).The density of the grating changes along theY-axis.The total number of gratings can be expressed as

    Fig.12.Diagram of the curved surface variable line pitch grating.

    where,N(y) is the total number of score lines accumulated from the originOto any pointPin theYdirection,d0is the score line spacing at the originO, and ζiis the spatial variation parameter of the score line density.Bringing the coordinates of each point into the optical path function ofAPBand expanding it into the Taylor series ofyandz, we obtain

    In Eq.(11),jandkcannot simultaneously be zero,and each term corresponds to a grating aberration of the curved grating.By optimizing the spatial variation parameters of the curved grating, the spectral focus curve and the spatial focus curve can be made to intersect at the astigmatic point, thereby achieving non-atigmatism imaging.Through the optimization of superposition polynomials,high-order aberrations can be corrected and the optimal solution of the system can be obtained.

    5.6.Computational Imaging Method

    The stitched image will inevitably have distortion,caused by the slit-bar data recorded with the CMOS image sensor, like a smile curve.Although small, this distortion is still non-negligible.This is caused by the inherent characteristics of the curved grating and off-axis imaging.Therefore, it is necessary to perform computational post-processing on the distorted stitched solar coronal image to obtain a corrected solar coronal image.The essence of the image distortion correction problem is to spatially transform the distorted image, and the transformation relationship is fixed in this system.A DL CNN can be used to fit this transformation[32-34], which is shown in Fig.13.The entire process can be viewed as feature extraction of the input image and pixel-by-pixel regression of the output results[35,36].Therefore, the CNN first convolves and down-samples the input data, performs feature transformation in 9 residual blocks, and finally up-samples through transposed convolution to reconstruct a distortion-free image.When training the model, we use the mean square error loss as the loss function, which is expressed as

    Fig.13.Schematic diagram of our deep learning CNN.The input is a distorted stitched solar coronal image and the output is a corrected solar coronal image.The network layers used in the CNN structure are shown to the right.

    where,Yiis the ideal coronal image,is the reconstructed coronal image, andnis the number of elements of the output tensor.

    6.AVAILABILITY OF DATA AND MATERIALS

    All data needed to evaluate the conclusions in this study are presented in this manuscript.Additional data related to this paper may be requested from the authors.

    7.ABBREVIATIONS

    ASMA: Annular Spherical Mirror Aperture

    CME: Coronal Mass Ejection

    CMOS: Complementary Metal-Oxide-Semiconductor

    CNN: Convolutional Neural Network

    COA: Circular Obscuration Aperture

    DL: Deep Learning

    EO: External Occulter

    EUV: Extreme Ultraviolet

    FOV: Field Of View

    MTF: Modulation Transfer Function

    RMS: Root Mean Square

    SDO: Solar Dynamics Observatory

    ACKNOWLEDEGMENTS

    This study is partially supported by the National Natural Science Foundation of China (NSFC) (62005120;62125504).

    AUTHOR CONTRIBUTIONS

    Cuifang Kuang and Xiangqun Cui conceived the idea and initiated the project.Jinming Gao, Yinxu Bian and Yue Sun mainly wrote the manuscript and produced the figures.Jinming Gao and Yinxu Bian conducted the optical design and data experiments.Jilong Peng and Qian Yu provided algorithm support.Xiangzhao Wang and Xu Liu edited the manuscript.Xu Liu and Xiangqun Cui supervised the project.All authors read and approved the final manuscript.

    DECLARATION OF INTERESTS

    Xiangqun Cui is the editor-in-chief for Astronomical Techniques and Instruments and was not involved in the editorial review or the decision to publish this article.The authors declare no competing interests.

    亚洲自拍偷在线| 欧美在线一区亚洲| 一级毛片女人18水好多| 日韩欧美三级三区| 国产激情久久老熟女| 无限看片的www在线观看| 欧美中文日本在线观看视频| 久久人妻福利社区极品人妻图片| 精品电影一区二区在线| 免费在线观看成人毛片| 97碰自拍视频| 在线观看午夜福利视频| 成人国产综合亚洲| 午夜精品在线福利| 亚洲精品在线观看二区| 久久久精品国产亚洲av高清涩受| 中文字幕人成人乱码亚洲影| 国内精品久久久久精免费| 天天躁狠狠躁夜夜躁狠狠躁| 欧美日韩亚洲国产一区二区在线观看| 亚洲熟妇中文字幕五十中出| 亚洲专区字幕在线| 免费在线观看视频国产中文字幕亚洲| 天天躁夜夜躁狠狠躁躁| 亚洲熟妇中文字幕五十中出| 国产伦一二天堂av在线观看| 1024香蕉在线观看| 啦啦啦观看免费观看视频高清| 一夜夜www| 国产精品98久久久久久宅男小说| 久久这里只有精品19| 久久人人精品亚洲av| 亚洲av电影不卡..在线观看| 精品久久久久久,| 欧美日韩亚洲国产一区二区在线观看| 他把我摸到了高潮在线观看| 热99re8久久精品国产| 一本精品99久久精品77| 男女做爰动态图高潮gif福利片| 国产成人av教育| 搞女人的毛片| 国产精品 国内视频| 变态另类成人亚洲欧美熟女| 欧美成人性av电影在线观看| 夜夜看夜夜爽夜夜摸| 成年免费大片在线观看| 久久人人精品亚洲av| 国内精品久久久久精免费| 国产亚洲欧美在线一区二区| 好看av亚洲va欧美ⅴa在| cao死你这个sao货| netflix在线观看网站| 亚洲五月婷婷丁香| 亚洲一区二区三区色噜噜| 色av中文字幕| 亚洲国产精品sss在线观看| 久久热在线av| 999精品在线视频| www国产在线视频色| 国内精品久久久久精免费| 少妇被粗大的猛进出69影院| 精品午夜福利视频在线观看一区| 国产av一区二区精品久久| 国产97色在线日韩免费| 搡老岳熟女国产| www.精华液| 久久狼人影院| 婷婷丁香在线五月| 一区福利在线观看| www.精华液| 美女高潮喷水抽搐中文字幕| 一本综合久久免费| 女警被强在线播放| 国产成人欧美| 人人妻人人澡人人看| 亚洲三区欧美一区| 国内精品久久久久久久电影| 国产私拍福利视频在线观看| 一级毛片女人18水好多| 99在线人妻在线中文字幕| 免费在线观看视频国产中文字幕亚洲| 久久亚洲真实| 悠悠久久av| 免费观看人在逋| 精品少妇一区二区三区视频日本电影| 国产成人av教育| 中文字幕另类日韩欧美亚洲嫩草| 亚洲国产精品sss在线观看| 黑人操中国人逼视频| 不卡一级毛片| 国产精品 欧美亚洲| 香蕉国产在线看| 老汉色av国产亚洲站长工具| 手机成人av网站| 精品日产1卡2卡| 欧美在线一区亚洲| 免费在线观看完整版高清| 黑人巨大精品欧美一区二区mp4| 欧美丝袜亚洲另类 | 中文字幕另类日韩欧美亚洲嫩草| 麻豆久久精品国产亚洲av| 黄色毛片三级朝国网站| 老司机午夜十八禁免费视频| 90打野战视频偷拍视频| 久久婷婷成人综合色麻豆| 天堂影院成人在线观看| 国产精品久久久久久精品电影 | 欧美又色又爽又黄视频| 人人澡人人妻人| 欧美人与性动交α欧美精品济南到| 免费电影在线观看免费观看| 黄色成人免费大全| 国产亚洲av嫩草精品影院| 天天躁狠狠躁夜夜躁狠狠躁| 免费一级毛片在线播放高清视频| 精品一区二区三区四区五区乱码| 香蕉丝袜av| 久久久久久九九精品二区国产 | 免费高清在线观看日韩| 超碰成人久久| 免费女性裸体啪啪无遮挡网站| 91字幕亚洲| 免费观看精品视频网站| 久9热在线精品视频| 久久亚洲精品不卡| 精品国产超薄肉色丝袜足j| 日韩欧美国产在线观看| 岛国在线观看网站| 99热这里只有精品一区 | 午夜久久久久精精品| 91av网站免费观看| 神马国产精品三级电影在线观看 | 亚洲午夜精品一区,二区,三区| 黄色a级毛片大全视频| 成年人黄色毛片网站| 亚洲国产欧洲综合997久久, | 欧美一级毛片孕妇| av视频在线观看入口| 精品不卡国产一区二区三区| 正在播放国产对白刺激| 大香蕉久久成人网| 国产黄色小视频在线观看| 免费在线观看完整版高清| 精品国产超薄肉色丝袜足j| 精品国内亚洲2022精品成人| 一级作爱视频免费观看| 亚洲国产欧洲综合997久久, | 12—13女人毛片做爰片一| 亚洲欧洲精品一区二区精品久久久| 欧美丝袜亚洲另类 | 一本一本综合久久| 一区二区日韩欧美中文字幕| 在线观看免费午夜福利视频| 国产av一区在线观看免费| 黑丝袜美女国产一区| 欧美日韩精品网址| 在线播放国产精品三级| 国产一区二区三区在线臀色熟女| 免费在线观看完整版高清| 老熟妇仑乱视频hdxx| 亚洲专区字幕在线| 色综合欧美亚洲国产小说| 午夜福利视频1000在线观看| 91成年电影在线观看| 91老司机精品| 久99久视频精品免费| 一区二区日韩欧美中文字幕| 亚洲 欧美一区二区三区| 午夜老司机福利片| 久久天躁狠狠躁夜夜2o2o| 中文字幕人妻丝袜一区二区| 日韩 欧美 亚洲 中文字幕| 999精品在线视频| 黄色成人免费大全| 久久国产精品人妻蜜桃| 大型av网站在线播放| 国产成人系列免费观看| 亚洲一区二区三区色噜噜| 亚洲午夜理论影院| 在线播放国产精品三级| 国产麻豆成人av免费视频| 香蕉丝袜av| 男女床上黄色一级片免费看| 97人妻精品一区二区三区麻豆 | 美女国产高潮福利片在线看| 日日爽夜夜爽网站| 亚洲人成网站高清观看| xxxwww97欧美| 99国产精品一区二区蜜桃av| 亚洲真实伦在线观看| 99国产精品99久久久久| 两个人看的免费小视频| 免费在线观看黄色视频的| 久久久水蜜桃国产精品网| 亚洲国产精品合色在线| 黄片播放在线免费| 十分钟在线观看高清视频www| 久久精品成人免费网站| 久久天堂一区二区三区四区| 国产高清视频在线播放一区| 欧美大码av| 日韩有码中文字幕| 亚洲精品国产精品久久久不卡| 日韩精品中文字幕看吧| 久久久久久久久免费视频了| 精品免费久久久久久久清纯| 成人特级黄色片久久久久久久| 国产熟女xx| 99热6这里只有精品| 人成视频在线观看免费观看| 成在线人永久免费视频| 国产野战对白在线观看| 黄色成人免费大全| 国产国语露脸激情在线看| 成在线人永久免费视频| 1024香蕉在线观看| 一进一出好大好爽视频| 国产99白浆流出| а√天堂www在线а√下载| 亚洲全国av大片| 精品久久久久久久久久久久久 | 每晚都被弄得嗷嗷叫到高潮| 高潮久久久久久久久久久不卡| 一个人观看的视频www高清免费观看 | 免费在线观看亚洲国产| 国产精品久久久av美女十八| 国产成人欧美| 国产视频一区二区在线看| 一级片免费观看大全| 久久精品aⅴ一区二区三区四区| 国产亚洲精品第一综合不卡| 国产伦在线观看视频一区| 两个人看的免费小视频| 国产久久久一区二区三区| 欧美乱妇无乱码| 日韩高清综合在线| 亚洲人成77777在线视频| 亚洲精品中文字幕在线视频| 欧美国产精品va在线观看不卡| 亚洲国产日韩欧美精品在线观看 | 色尼玛亚洲综合影院| 欧美色视频一区免费| 老汉色av国产亚洲站长工具| 国产99白浆流出| 免费无遮挡裸体视频| 欧美大码av| 国产日本99.免费观看| 亚洲五月天丁香| 一区二区日韩欧美中文字幕| 久久久久免费精品人妻一区二区 | 伦理电影免费视频| av欧美777| av免费在线观看网站| 午夜成年电影在线免费观看| 亚洲成av人片免费观看| 一夜夜www| 亚洲成av片中文字幕在线观看| 欧美日本亚洲视频在线播放| 亚洲精品中文字幕在线视频| 国产激情欧美一区二区| 欧美中文日本在线观看视频| a在线观看视频网站| 国产精品久久久久久亚洲av鲁大| 中文字幕高清在线视频| 一级作爱视频免费观看| 欧美丝袜亚洲另类 | 国内精品久久久久久久电影| 欧美性猛交╳xxx乱大交人| 88av欧美| 两个人免费观看高清视频| 美女午夜性视频免费| 日韩欧美国产一区二区入口| 老司机靠b影院| 女性被躁到高潮视频| 搞女人的毛片| 中文资源天堂在线| av片东京热男人的天堂| 午夜久久久在线观看| 精品一区二区三区四区五区乱码| 国产主播在线观看一区二区| 成年版毛片免费区| 国产成人一区二区三区免费视频网站| а√天堂www在线а√下载| av在线天堂中文字幕| 日韩三级视频一区二区三区| 国产爱豆传媒在线观看 | 国产又黄又爽又无遮挡在线| 人妻久久中文字幕网| 天堂影院成人在线观看| 国内精品久久久久精免费| 亚洲国产精品久久男人天堂| 老司机午夜福利在线观看视频| 亚洲欧美精品综合久久99| 亚洲精品国产一区二区精华液| 国产亚洲欧美精品永久| 美女扒开内裤让男人捅视频| 国产亚洲欧美精品永久| 亚洲av电影不卡..在线观看| 999久久久国产精品视频| 999久久久国产精品视频| 午夜福利18| 一区二区日韩欧美中文字幕| 久久狼人影院| 男女做爰动态图高潮gif福利片| 国产精品香港三级国产av潘金莲| 亚洲久久久国产精品| 日韩三级视频一区二区三区| 免费观看人在逋| 51午夜福利影视在线观看| 视频区欧美日本亚洲| 国产精品久久久人人做人人爽| 97人妻精品一区二区三区麻豆 | 亚洲一区二区三区不卡视频| 怎么达到女性高潮| 欧美精品啪啪一区二区三区| 精品福利观看| 亚洲av五月六月丁香网| 午夜福利高清视频| 91字幕亚洲| 在线视频色国产色| 男人舔女人下体高潮全视频| 男女做爰动态图高潮gif福利片| 久久精品国产清高在天天线| 亚洲精品美女久久久久99蜜臀| 亚洲精品美女久久久久99蜜臀| 国产午夜精品久久久久久| 国产一区二区在线av高清观看| 满18在线观看网站| 久久久久久久精品吃奶| 中文字幕另类日韩欧美亚洲嫩草| 大型av网站在线播放| 国产精品av久久久久免费| 日韩欧美免费精品| 日本免费a在线| 欧美性猛交╳xxx乱大交人| 日韩高清综合在线| 麻豆成人午夜福利视频| 国产亚洲欧美98| 法律面前人人平等表现在哪些方面| 91老司机精品| 一级毛片女人18水好多| 久久久久国产精品人妻aⅴ院| 在线观看免费视频日本深夜| 一区二区三区高清视频在线| 波多野结衣高清无吗| 日韩欧美国产一区二区入口| 亚洲,欧美精品.| 久久久久久大精品| av视频在线观看入口| a在线观看视频网站| 一区二区三区精品91| 别揉我奶头~嗯~啊~动态视频| bbb黄色大片| 国产成人一区二区三区免费视频网站| 成人18禁在线播放| 久久精品国产综合久久久| av免费在线观看网站| 妹子高潮喷水视频| 久久久久久人人人人人| 欧美激情 高清一区二区三区| 国产欧美日韩精品亚洲av| 欧美成人性av电影在线观看| 国产精品爽爽va在线观看网站 | 成人av一区二区三区在线看| 亚洲一区二区三区色噜噜| 伊人久久大香线蕉亚洲五| 12—13女人毛片做爰片一| 午夜视频精品福利| 成人国产一区最新在线观看| 少妇粗大呻吟视频| 精品国产超薄肉色丝袜足j| 亚洲 国产 在线| 男女午夜视频在线观看| 欧美av亚洲av综合av国产av| avwww免费| 亚洲片人在线观看| 视频区欧美日本亚洲| 搡老妇女老女人老熟妇| 国产伦人伦偷精品视频| 在线观看午夜福利视频| 一本一本综合久久| 女人被狂操c到高潮| 欧美中文综合在线视频| 欧美最黄视频在线播放免费| 超碰成人久久| 亚洲国产欧美日韩在线播放| 男人舔奶头视频| 在线观看66精品国产| 亚洲专区国产一区二区| 每晚都被弄得嗷嗷叫到高潮| 99re在线观看精品视频| av欧美777| 久久婷婷人人爽人人干人人爱| 国产精品九九99| 俄罗斯特黄特色一大片| 91国产中文字幕| 怎么达到女性高潮| 国内少妇人妻偷人精品xxx网站 | 精品久久久久久久久久久久久 | 高清在线国产一区| 欧美亚洲日本最大视频资源| 国产精品精品国产色婷婷| 欧美丝袜亚洲另类 | 国产伦人伦偷精品视频| 满18在线观看网站| 又大又爽又粗| 国产日本99.免费观看| 十分钟在线观看高清视频www| 亚洲国产精品久久男人天堂| 国产成人一区二区三区免费视频网站| 中文字幕人妻熟女乱码| 久久草成人影院| 中文字幕最新亚洲高清| 国产激情欧美一区二区| 免费在线观看视频国产中文字幕亚洲| 午夜两性在线视频| 国产男靠女视频免费网站| 国产一区二区三区在线臀色熟女| 黄片大片在线免费观看| 成年女人毛片免费观看观看9| 国产片内射在线| 我的亚洲天堂| 久久香蕉精品热| 成人国产一区最新在线观看| av欧美777| 中文字幕人妻丝袜一区二区| 亚洲熟妇中文字幕五十中出| 国产欧美日韩精品亚洲av| 91国产中文字幕| 欧美成人免费av一区二区三区| 久久热在线av| 精品卡一卡二卡四卡免费| 国产精品国产高清国产av| avwww免费| 黄色视频不卡| 欧美成人一区二区免费高清观看 | 女警被强在线播放| 欧美黄色淫秽网站| 国产成+人综合+亚洲专区| 欧美日韩黄片免| 精品国产一区二区三区四区第35| 亚洲国产中文字幕在线视频| 午夜精品久久久久久毛片777| 一本精品99久久精品77| 亚洲av片天天在线观看| 搞女人的毛片| 久久久久久人人人人人| 日韩精品青青久久久久久| av在线播放免费不卡| 淫秽高清视频在线观看| 后天国语完整版免费观看| 一本精品99久久精品77| 欧美不卡视频在线免费观看 | 精品卡一卡二卡四卡免费| 国产精品一区二区精品视频观看| 黄片小视频在线播放| 国产精品久久视频播放| 中文在线观看免费www的网站 | 免费在线观看成人毛片| 欧美另类亚洲清纯唯美| 不卡一级毛片| 欧美黄色片欧美黄色片| 亚洲精品中文字幕一二三四区| 国产97色在线日韩免费| 91九色精品人成在线观看| 黑人巨大精品欧美一区二区mp4| 色综合亚洲欧美另类图片| 免费看十八禁软件| 国产单亲对白刺激| netflix在线观看网站| 色哟哟哟哟哟哟| 19禁男女啪啪无遮挡网站| 精品无人区乱码1区二区| 亚洲天堂国产精品一区在线| 99久久99久久久精品蜜桃| 熟女电影av网| videosex国产| 日韩高清综合在线| 天堂影院成人在线观看| 日韩欧美 国产精品| 国产精品日韩av在线免费观看| 久久精品成人免费网站| 女人被狂操c到高潮| 日本成人三级电影网站| 丰满的人妻完整版| 99久久无色码亚洲精品果冻| 国产视频内射| 黑人巨大精品欧美一区二区mp4| 久久亚洲精品不卡| 欧美激情高清一区二区三区| 1024手机看黄色片| www国产在线视频色| 在线播放国产精品三级| 婷婷亚洲欧美| 国产欧美日韩精品亚洲av| 黄片大片在线免费观看| 日韩国内少妇激情av| 免费看a级黄色片| 听说在线观看完整版免费高清| 欧美另类亚洲清纯唯美| 女性被躁到高潮视频| 亚洲成人久久性| 亚洲七黄色美女视频| 国产v大片淫在线免费观看| 黄网站色视频无遮挡免费观看| 欧美性猛交黑人性爽| 两人在一起打扑克的视频| 国语自产精品视频在线第100页| 哪里可以看免费的av片| 一本一本综合久久| 久久青草综合色| 国产又爽黄色视频| av超薄肉色丝袜交足视频| 成熟少妇高潮喷水视频| 国产精品免费视频内射| 91成年电影在线观看| 在线看三级毛片| 在线观看舔阴道视频| 亚洲第一av免费看| 亚洲国产欧美一区二区综合| 可以免费在线观看a视频的电影网站| 欧美乱色亚洲激情| 亚洲av成人不卡在线观看播放网| 久久久久久久午夜电影| 国产成人一区二区三区免费视频网站| 亚洲午夜理论影院| 村上凉子中文字幕在线| 亚洲五月天丁香| 欧美国产日韩亚洲一区| 在线观看66精品国产| 久久亚洲精品不卡| 99riav亚洲国产免费| 久久热在线av| 妹子高潮喷水视频| 十八禁网站免费在线| 国产av又大| 免费在线观看日本一区| 男人舔奶头视频| 久久久国产精品麻豆| 可以在线观看的亚洲视频| 一级毛片高清免费大全| 午夜免费成人在线视频| 国产黄a三级三级三级人| 97超级碰碰碰精品色视频在线观看| 色综合婷婷激情| 国内少妇人妻偷人精品xxx网站 | 国产成人系列免费观看| 极品教师在线免费播放| 精品久久蜜臀av无| 亚洲专区字幕在线| 国内少妇人妻偷人精品xxx网站 | 久久久久久国产a免费观看| 制服丝袜大香蕉在线| 十分钟在线观看高清视频www| 国产爱豆传媒在线观看 | 亚洲精品国产区一区二| 免费电影在线观看免费观看| 中文字幕人妻熟女乱码| 精品高清国产在线一区| 欧美午夜高清在线| 桃色一区二区三区在线观看| 亚洲国产高清在线一区二区三 | 午夜影院日韩av| 啪啪无遮挡十八禁网站| www.精华液| 国产精品99久久99久久久不卡| 欧美激情极品国产一区二区三区| 国产单亲对白刺激| 欧美乱码精品一区二区三区| 亚洲第一av免费看| 两个人视频免费观看高清| 亚洲avbb在线观看| 亚洲 欧美 日韩 在线 免费| 国产午夜精品久久久久久| 热99re8久久精品国产| 999久久久国产精品视频| 国内毛片毛片毛片毛片毛片| 精品久久久久久久人妻蜜臀av| 一进一出抽搐gif免费好疼| 精品国产美女av久久久久小说| 久久久国产成人免费| 热re99久久国产66热| 亚洲精品美女久久av网站| 日韩一卡2卡3卡4卡2021年| 国产熟女xx| 97碰自拍视频| 香蕉国产在线看| 午夜日韩欧美国产| 欧美日韩亚洲综合一区二区三区_| 99久久99久久久精品蜜桃| 日韩视频一区二区在线观看| 欧美日本视频| 亚洲狠狠婷婷综合久久图片| 一区二区三区国产精品乱码| 免费看美女性在线毛片视频| 国产成+人综合+亚洲专区| 亚洲七黄色美女视频| 手机成人av网站| 高清毛片免费观看视频网站| 精品一区二区三区视频在线观看免费| 给我免费播放毛片高清在线观看| 妹子高潮喷水视频| 人人妻人人澡人人看| 人人妻,人人澡人人爽秒播| www国产在线视频色| 男女做爰动态图高潮gif福利片| 91大片在线观看| 波多野结衣高清作品| 中文字幕高清在线视频| 一个人免费在线观看的高清视频| 又黄又粗又硬又大视频| 欧美成人性av电影在线观看| 婷婷亚洲欧美| 啦啦啦 在线观看视频| 欧美久久黑人一区二区|