• 美文
  • 文章
  • 散文
  • 日记
  • 诗歌
  • 小说
  • 故事
  • 句子
  • 作文
  • 签名
  • 祝福语
  • 情书
  • 范文
  • 读后感
  • 文学百科
  • 当前位置: 柠檬阅读网 > 范文 > 正文

    高分辨率井周储层深探测反射声波测井成像方法

    时间:2023-03-10 08:15:05 来源:柠檬阅读网 本文已影响 柠檬阅读网手机站

    岳文正,田 斌

    (1.油气资源与探测国家重点实验室,中国石油大学(北京),北京102249;2.中国石油天然气股份有限公司辽河油田分公司沈阳采油厂测试公司,辽宁新民110316)

    随着油气勘探的不断深入和发展,复杂储层已成为研究的重点。在复杂储层评价中,测井和地震勘探是两种关键的地球物理方法,综合利用测井资料和地震资料可以获取丰富的地层信息。但是,由于两种方法的物理机理差异,它们的分辨率和探测范围明显不同,故而很难直接将两种方法的测量结果进行匹配来识别储层。因此,研究能够把测井和地震数据结合起来的成像方法是解决该问题的关键之一。目前,基于阵列声波测井的反射波成像方法是实现这一目标的主要发展方向。

    反射波成像主要利用辐射到井外的能量,对井周地层进行成像[1-4]。利用信号处理方法提取反射波,结合偏移成像技术进行井旁构造成像。自20世纪80年代以来,该方法在国内外测井领域快速发展[3-8],并且在裂缝性[7]以及缝洞型[8]碳酸盐岩储集层评价中显示出突出的应用效果。其中,以探测井周反射体方位和距离为主要目标的远探测反射声波成像研究最为广泛。

    远探测声波测井相关研究起步较早,速度与衰减评价仪器(EVA)是最早的远探测声波测井原型机[9]。1989年,HORNBY[1]将反射波从测井全波列数据中提取出来,实现了井周探测成像,标志着声波远探测技术的诞生。1998年,在偶极横波测井仪(DSI)基础上研发的反射波成像测井仪(BARS)开启了声波远探测技术的商业化,而2006年研制的Sonic Scanner仪器进一步推广了远探测声波测井的应用[2,10]。此外,针对多极子阵列声波测井仪器(XMAC Ⅱ)的声场特点,开发了相应的成像处理算法[11-12]。2004年,TANG[13]利用四分量交叉偶极数据对井下构造进行成像,发现低频偶极声源产生的偶极反射波比传统单极反射波具有更大的探测深度,且具有一定的方位识别能力。

    国内相关领域的研究始于理论分析[14-21]。1998年,王乃星等[14]针对声波测井的全波列信号研究了反射波的分离与成像处理,并计算反射界面视倾角的影响。2001年,车小花等[15]进行了模拟实验,以水和PVC板的交界面模拟裂缝,使用中值滤波进行反射波的分离,根据反射波列计算裂缝倾角。2002年,薛梅等[16]从几何声学出发,探究多种因素对反射波和滑行波接收结果的影响,为反射声波测井仪器的研制提供了理论基础。2003年,李长文等[17]分析了不同裂缝发育状态的测井响应,据此开发针对反射斯通利波、反射纵波的处理解释方法,并探讨了反射波信息的进一步应用。2004年,宋立军等[18]建立柱状三层介质井孔体系物理模型,实现了超长源距声波测井全波数值模拟与分析,并得到了远井界面时间记录剖面成像图。2005年,何峰江等[19]利用三维有限差分数值模拟方法研究了声反射成像测井的波场特征,基于小波变换来压制噪声和提高反射波信噪比,并在反射波信号处理中引入了多尺度相关法。2012年,WEI等[20]和唐晓明等[21]研究了偶极子声源的远场辐射特征,开发偶极横波成像方法,认为SH横波是偶极横波远探测成像的最佳信号。除算法外,针对传统反射声波测井技术周向分辨能力不足的问题,乔文孝等[22]研究了声波换能器的性能,设计了相控声波技术,可以实现三维反射声波的成像。

    远探测声反射成像测井在很多领域取得了明显的应用效果,这种方法主要探测井周反射体的位置和方位。但是,为了更直观地显示井周储层的岩性、流体、微构造和裂缝、溶洞等,需要开发新的成像方法。近年来,一种基于高分辨率成像的深探测反射波成像方法逐渐发展起来[23-24],该方法主要对井周岩性、流体分布、微构造变化等进行成像,由于其探测深度为数十米,分辨率为十几厘米,所以初步实现了测井-地震联合成像储层识别分析。

    本文介绍了高分辨率深探测反射波成像新方法及其在井周复杂储层反射体、流体分布评价等领域的应用。该方法利用Shearlet变换提取阵列声波测井反射波,基于层析速度成像获取井周径向速度剖面,在此基础上结合偏移算法对井周地层进行高分辨率多信息成像。利用该方法对我国某油气区井数据进行处理,开展测井-地震联合成像分析,对井周岩性、流体分布等进行成像,并分析压裂前后地层的改造情况。

    图1展示了声波反射成像测井的基本原理,基于传统或改进的阵列声波仪器接收井下地层中的弹性波,其中包括井壁直达波(滑行波)和远井眼地层反射波。

    图1 声波反射成像测井原理示意

    高分辨率深探测声波测井成像,主要基于阵列声波测井数据,利用井中反射波对井周岩性、流体分布等信息进行成像处理,所以,如果反射波提取效果不理想,那么成像结果在显示井周反射体的细节信息方面会存在不足。因此,该方法的关键在于反射波的提取。

    1.1 反射波提取方法

    阵列声波测井数据的主要特点是:①仪器接收到的反射波信号幅度相比模式波而言微弱得多,并且在某些源距下,反射波和模式波会相互混叠,常规的反射波提取方法很难实现波的有效分离;②对反射信号的处理,需要做上、下行反射波的分离。因此,高效和适用性强的反射波提取方法一直是众多学者研究的热点。

    F-K滤波法最早被用于从全波列数据中提取反射信号[1]。LI等[11]基于VSP处理相关方法,在波场分离中综合运用了中值滤波和F-K滤波。TANG[13]为了解决单极声源激发无方位识别能力的问题,提出了一种指向性测量方法进行远探测成像,采用偶极激发和偶极接收,综合使用了F-K滤波和中值滤波,以四分量数据的反射波进行成像,获得了构造的方位信息。TAO等[25]提出了多尺度相关法提取反射信号,主要结合双树复小波(DT-CWT)变换与慢度-时间相关(STC)方法分离出反射波。范宜仁等[26]针对XMAC实测数据,运用多尺度相关法成功分离出了反射波。王兵等[27]将多尺度相关法与参数估计法[28]和F-K滤波进行比较,通过对比模拟数据和实际数据的应用效果,证实了多尺度相关法的优越性。但是,多尺度相关法存在的问题是参数控制较为麻烦,无法对整个地层统一处理,某分块处理的方式影响了该方法的效率和推广应用[29]。刘喜武等[30]和宋建国等[31]推导了基于求解高分辨率Radon变换的稀疏约束共轭梯度算法,并将其用于多次波的压制。李超等[32]将这一方法引入到反射波的提取中,以阻尼最小二乘解作为初始解进行高分辨率Radon变换的迭代求解,实现了反射波提取。实际数据处理结果表明,高分辨率Radon变换在反射波分离中的应用效果良好[33]。

    针对反射波提取,虽然算法较多,但是由于井孔中接收到的信号其滑行波能量强,反射波信号微弱,且具有多种模式[19],故常淹没在滑行波中难以分离,所以大多数的方法未能完全解决波的分离问题。非线性滤波方法[34],具有保留指定方向主要数据的特性,前人的研究[11,35]主要是针对某一道集进行处理,并没有指定适用条件,且单一的滤波难以适应噪声严重的实际数据。针对这些问题,2019年,CHEN等[36]将Shearlet变换引入到阵列声波测井的波场分离中,选取共接收道集进行Shearlet变换处理。

    Shearlet变换[37-39]是一种近似最优函数稀疏表示方法,不同特征的Shearlet函数则是通过对基函数进行缩放、剪切及平移等仿射变换处理得到。Shearlet变换与小波变换具有相似性,是小波变换的进一步优化后的折变换法。引入Shearlet变换,实现了信号的更高精度多尺度分解。此外,Shearlet变换具有灵活的方向选择性,能对信号进行稀疏表示并产生最优逼近。

    深探测反射波提取是基于Shearlet变换利用Meyer小波函数和冲击函数计算母函数,再通过母函数的膨胀、方位剪切以及时间域平移得到基函数,利用信号和基函数的内积,可以计算得到不同尺度、方位以及平移的变换结果,将信号转换到剪切波域。定义双曲尺度矩阵Mp和剪切矩阵Ms为:

    (1)

    (2)

    式中:p为膨胀参数;s为剪切参数。Shearlet变换的基函数可以定义为:

    (3)

    式中:b为平移参数。如果将基函数从时间域变换到傅里叶域,就可以得到公式(4)。其中母函数是由两个一维函数之积来计算(公式(5)),函数ψ1表示Meyer小波,函数ψ2代表冲击函数。

    (4)

    (5)

    这样,信号f的Shearlet变换可以用公式(6)来表示。

    (6)

    基于Parseval格架理论,可以得到变换结果的逆变换,该逆变换可以利用正变换得到的剪切系数重构原信号[40]。观察阵列声波测井数据的共偏移距道集,可以发现不同波的到时特征。通常,模式波的到时表现为线性特征,而反射波表现为曲线特征。在Shearlet域,模式波能量将聚集在小尺度的水平扇形区域,因此只需利用阀值便可以将模式波分离出来。为了验证该方法的有效性,采用交错网格有限差分的正演模拟合成声波测井数据,考察不同反射波提取方法的应用效果,提取结果如图2所示。

    图2 各种反射波提取方法的提取结果a 原始波形数据; b 理论反射波形; c Shearlet变换提取结果; d Radon变换提取结果; e 小波变换提取结果; f 中值滤波提取结果

    从图2可以看出,相较于Radon变换、小波变换和中值滤波的提取结果,Shearlet变换提取的反射波形更加完整,能量损失少,剩余残差干扰小。为了进一步说明,选出不同方法提取结果中的同一道反射波数据(图3中标为红色),与理论反射波形进行对比。图3中的相关系数进一步表明,Shearlet变换方法提取反射波的效果明显好于其它方法,可以提取高质量反射波数据。

    图3 理论反射波形与提取的反射波对比(R为波形相关系数)a Shearlet变换提取结果; b Radon变换提取结果; c 小波变换提取结果; d 中值滤波提取结果

    1.2 井周速度层析成像

    在有效提取反射波信号的基础上,利用偏移算法对井周储层进行成像处理是深探测声反射成像的又一关键所在。对于偏移成像方法,影响其成像质量的最重要因素是速度模型。一般而言,可以认为某一深度井周地层的速度是均匀的,即井周的地层速度与远井眼地层速度相同。在很多研究中,地层的速度模型均假设为某一深度处一致。然而,实际钻井过程可能会产生裂缝和井周地层的泥浆侵入,这些因素会导致在井眼周围形成变化的速度剖面。为了得到高质量成像结果,在本研究中利用速度层析成像获取和分析井周的速度剖面,建立符合井周分布的速度模型,而不是采用常规的均一层状速度模型。

    在井周速度层析成像方面,早期的方法主要利用折射初至波,采用迭代重建技术进行层析成像[41-42],但是这种方法计算量大,难以有效应用。夏克文等[43]、李昌彪等[44]分别探讨了利用反射初至波进行层析成像的方法。WINKLER[45]和BALLAND等[46]

    探讨了利用超声波探测重建井周速度场的方法。SMAINE等[47]根据滑行波到时信息反演了井周二维速度剖面。TANG等[48]提出了基于偶极横波频散曲线的横波速度径向剖面反演方法。

    纵波信号从单极子声源T发出,经过速度不同的i个侵入层滑行传播到i个接收器中被接收到。从图4a 中可以看到,短源距能接收到浅地层的滑行纵波信息,而长源距可以接收到较远地层的滑行纵波信息。测井记录下滑行波到达每个接收器的时间,从这i个到时信息中可以反演出井壁附近侵入层的速度和厚度信息。如果每个深度点有M个不同源距的接收器,那么每个深度点就会对应M个源距接收到的M个不同到时信息,则假设每个深度点有M个速度层,利用这些不同接收器组合接收到的到时信息,结合公式(7)就可以反演出井壁附近M个不同的速度层。

    (7a)

    (7b)

    (7c)

    式中:vf为井孔内泥浆纵波速度;s为仪器源距;在井壁附近形成i层速度层数,每个层的厚度为Hi,相对应的速度表示为vi。图4给出了某油田井周不同方向的速度变化剖面。

    图4 深探测声波测井示意以及井周不同方向速度变化剖面a 测井示意; b 北向(N); c 东向(E); d 南向(S); e 西向(W)

    图4b至图4e给出了某油田井周不同方向的速度变化剖面,横坐标为距井眼距离。图中蓝色表示速度未发生变化区间,红色部分代表速度变化较大的区间。图中结果显示,井周不同深度井周速度变化范围不同。此外,井周不同方向速度变化亦不同,结合这些速度的变化特征,可以对地层的各向异性进行评价。

    1.3 高分辨率偏移成像

    偏移的目的是将倾斜反射归位到真正的位置,并使绕射波收敛,得到真实的地下构造图像。深探测反射声波测井成像的偏移成像方法一般借鉴较为成熟的地震资料处理方法,并且针对声波测井数据的特殊情况,选用较为合适的叠前偏移方法,再进行参数调整,得到最终的偏移成像结果。

    1989年,HORNBY[1]在反射波成像处理中使用了广义Radon变换,通过共中心点道集叠加压制噪声,以提高反射波质量,刻画井旁构造信息。其后,BORLAND等[49]、HALDORSEN等[50]和MAIA等[51]采用周向排列的阵列接收器,得到不同方位接收的声波阵列数据,然后使用广义Radon变换偏移方法对深度点所有方位都进行偏移成像,再进行数据的方位聚焦,得到相对应的裂缝和地层的信息。ESMERSOY等[2]使用叠前Kirchhoff偏移处理反射声波数据,计算出井眼轨迹与地层的相对位置关系,其计算结果与实际情况吻合较好。在此之后,COATES等[3]、LI等[11]和TANG等[6,33]都使用Kirchhoff积分偏移开展相应的成像处理研究,取得了很好的应用结果。

    大多数偏移方法的使用都基于分离后的上下行反射波场数据,而叠前F-K偏移的一大优势在于仅需进行一次偏移,并且偏移后直接分离出上下行波场的结果,同时该偏移方法没有倾角的限制。但是F-K偏移的基础是三维数据,因此针对测井数据的结构需要进行相应的改进。ZHENG等[12]根据反射声波测井的相关特点,提出了改进的F-K偏移方法。魏周拓等[52]通过建立的高角度井旁地质界面模型,基于模拟数据验证了改进的叠前F-K偏移成像方法的适用性,并针对XMACⅡ仪器的实测数据,证实了该方法可以有效确定井旁地层界面的相关信息。

    逆时偏移(reverse time migration,RTM)是近年来地面地震主流的偏移成像技术,对复杂情况下的成像处理效果较好。该方法基于完整的双程波波动方程,可以有效地对复杂情况下的接收波场进行有效成像。相较于传统成像方法,逆时偏移的成像结果更准确、信噪比更高、对于地质构造体的刻画也更加清晰。逆时偏移的首次提出时间比较早,但是由于其要求的计算量和存储量相当巨大,早期的应用受到了很多限制。而随着近年来计算机计算能力的不断提升,逆时偏移技术也得到极大的发展。众多学者对逆时偏移的成像条件、储存策略及压制噪声等方面开展了深入研究,但就其方法对深探测声波测井井孔环境的适应性上,还需要做进一步的完善[28,53-56]。

    2.1 测井-地震联合成像

    在我国某探区,其油气勘探的特点是发现断层即发现储层。因此,该区储层探测主要是以断层为主要勘探目标。图5给出该区某一探测结果,根据地震成像结果(图5a),可以判断有2个断层发育(红色曲线标示),为了评价断层周围储层发育情况,钻了大斜度井,井眼轨迹用绿色曲线在图5a中标出,可以看出X738m深度处,井周发育一孤立反射体。对该井进行了阵列声波测井,利用声波测井数据提取反射波进行高分辨率深探测成像处理,反射面结果如图5b所示,成像范围为图5a中矩形框区域Ⅰ。与地震成像对比可以发现,反射波测井成像结果显示孤立反射体的形状、位置、发育特征等与地震成像相同,而且,能够揭示异常体内部更小尺度的构造。因此,深探测反射波测井成像结果可以细化地震成像单一同相轴的内部小微构造信息,利用测井-地震联合成像可以给出更高分辨率的储层评价。除了反射面信息,深探测成像还可以给出岩性变化信息(图5c),可以看出图5c 中的岩性变化特征与图5a和图5b相一致。

    图5 同一数据的不同方法成像结果a 地震剖面; b 深探测反射成像; c 深探测岩性变化成像

    比较图5a和图5b发现,深探测方法具有较高的成像质量和分辨率,相较于常规方法在显示井周18m范围内反射体的细节信息方面具有优势。因此,这些结论表明深探测成像方法是一种有效的反射波提取和成像手段。结合地震剖面,可以分析反射体的细节信息,从而提高和精化复杂地层油气勘探的储层评价质量。

    2.2 页岩水平井地层岩性成像

    除了对反射体界面进行成像外,岩性和构造成像也是利用测井反射波深探测成像评价储层的重要研究方向。图6展示了一个对水平井井周构造进行深探测成像应用的实例,其中的测井数据源自我国西部的典型页岩气储层。图6a显示了利用地震成像解释获取的该地区基本地质构造情况,其中,井眼轨迹标注为红色光滑曲线。为了更清楚地显示结果,选取图6a 中Ⅰ和Ⅱ两个矩形框内的成像处理结果放大展示,分别如图6b和6c所示,图中不同颜色代表地层岩性。

    图6 水平井井周地层构造成像a 地震地质解释结果; b 图6a矩形框Ⅰ的深探测成像结果; c 图6a矩形框Ⅱ的深探测成像结果

    从图6中可以得出以下结论:①相较于图6a中的地震地质解释结果,深探测成像的结果可以提供井周20m范围内小尺度构造的更多细节信息,如图6b和图6c所示;②深探测成像方法的地层界面成像结果与地震地质解释结果具有一致性,可以直观观测和分析地层走向;③地层的局部不整合面和垮塌带可以从深探测成像结果中直接进行分析和解释,如图6b 中区域(1)的矩形框内数据所示;④在图6b区域(2)中,发育一条有限延伸的裂缝,并未完全穿层,因此没有导致明显的断层位移;⑤地震成像中断层的细节信息可以由深反射成像结果直观显示,如图6c中的区域(3)为断裂带,其对应于地震成像图6a中的断层区域Ⅱ,此外,图6c中的区域(1)和区域(2)显示了发育规模较小的两个断层,其在地震成像结果图6a 中未有显示。显然,所有结论表明,深探测成像方法能够实现井周反射体的高质量成像,并提供小尺度反射体的直观细节信息从而细化地震成像结果,可以用于页岩等复杂储层评价与识别。

    2.3 井周流体分布成像

    流体识别和量化一直都是地球物理储层评价的重点关注目标之一。深探测成像在提取反射界面信息的同时,进一步对井周的流体分布情况进行成像。该方法在国内多个油田以及海外土库曼斯坦等地区进行了现场试验和应用,现场应用识别气层有效性突出,计算的气体产能与单簇射孔结果吻合度高,最大应用井深超过8000m。图7为我国南部的某气田井周流体分布成像结果,其中图7a 中的气产能指数为基于深探测成像提取的气产能量化参数,数据源自于该气田某一页岩气层。图7b 和图7c中的数据源自该气田气层。

    图7 气产能指数及井周流体分布成像结果a 气产能指数对比验证; b 气产能指数曲线; c 井周流体分布成像结果

    图7a第4道中黄色充填曲线为气产能指数,红色杆状数据点为单簇射孔气测产能测试结果。分析图7a中的结果可知,深探测气产能指数与单簇射孔测试产能数据十分吻合,验证了产能指数的可靠性。图7c为井周流体分布成像结果,白色为气体。成像显示,气层顶界在X600m附近,主力产气层顶界在X650m,气油界面在X810m附近,于X880m处开始,产液量明显增大。深探测声波测井成像分析结果与测试结果吻合,进一步证实了该方法在储层含气层评价方面的适用性。

    2.4 地层压裂效果成像

    在页岩气等非常规油气勘探开发中,储层改造技术(包括压裂、酸化等)是提高储层产能的主要方法。但是,如何准确评价储层改造的效果一直缺少一种直观的手段。利用深探测成像的高分辨率能力,初步对储层压裂效果开展评价研究。图8显示了我国西北某气田一口井的深探测成像结果,该井分别在裸眼井阶段、套管井阶段、压裂完井阶段进行阵列声波测井,基于不同阶段测井数据,利用深探测成像方法获取井周高分辨率成像,直观评价储层改造效果。

    图8 同一口井不同阶段深探测成像结果a 裸眼井阶段; b 套管井阶段; c 压裂后阶段

    图8a至图8c分别为裸眼井、套管井和压裂后井周高分辨率成像结果。图中显示裸眼井成像结果和套管井成像结果保持一致,这说明:一方面是该方法具有稳定性,如图中的区域(3)在图8a至图8c中均未发生明显改变;另一方面是套管和水泥环对深探测成像的影响不明显。从图8c可以看到,压裂过程对储层产生了影响,尤其是区域(1)和区域(2),压裂后的储层地层特性发生明显变化。仔细对比图8a和图8c 可以直观评价压裂对储层的改造波及范围和程度,进而分析压裂的效果。因此,基于阵列声波测井的深探测成像方法为直观评价储层压裂改造效果提供了一种有效的手段。

    本文介绍了一种深探测反射声波测井成像技术及其在储层识别中的应用。该方法利用反射波而不是直达波(折射波)信息,能够对井旁裂缝、地质构造进行成像,为综合地质研究和油气勘探开发提供技术支撑。该方法在准确提取高分辨率反射波的基础上,通过井周速度层析成像和优化偏移算法,有效压制伪界面的影响,而且极大地改善了垂直于井眼的反射界面成像效果。利用该方法对实际阵列声波测井资料进行成像处理,建立断层、裂缝等反射体以及流体分布、岩性变化等识别特征,对井周的地质构造进行评价,进而识别储层。结合前文的讨论和分析,虽然该方法具有良好的应用效果,但是仍然需要在以下几个方面开展进一步的研究和探讨。

    1) 反射方位信息的获取和成像。深探测反射波成像虽然可以提供高分辨率的多信息成像,但是其结果多为二维的截面成像,无法提供截面延伸的方位信息。在油气资源勘探开发评价中,掌握储层的缝、洞以及流体分布等发育方位具有重要意义,尤其对于酸化压裂等增产措施,可以有效指导压裂方位的确定。

    2) 高分辨率井周三维成像。测井仪器和数据采集技术的发展,已经使得获取井周不同方向传输的反射波信号成为可能,因此如何充分利用不同方向反射波的丰富信息,对井周储层进行高分辨率的三维成像是将来研究的重要方向之一。如果结合深探测反射波成像方法,对井周反射体、岩性和流体分布进行三维分析评价,将更有利于复杂储层的精细描述,为深地油气资源探测和开发提供有力的新工具。

    3) 沿井壁传播反射波的自适应滤除。由于沿井壁传播的波,在层界面、过井裂缝等位置产生沿井壁传播的反射波,该信号并未携带地层深处的有用信息,并且随地层速度的不同而发生变化,对反射波成像产生干扰,因此,开发基于不同地层波速的自适应井壁反射波滤除方法也是需要研究的问题。

    猜你喜欢 高分辨率测井储层 本期广告索引测井技术(2022年3期)2022-11-25——北美又一种非常规储层类型">输导层
    ——北美又一种非常规储层类型石油与天然气地质(2021年5期)2021-10-29基于边界探测技术的储层识别与追踪化工管理(2021年7期)2021-05-13高分辨率合成孔径雷达图像解译系统雷达学报(2020年3期)2020-07-13基于储层构型研究的储层平面非均质性表征西南石油大学学报(自然科学版)(2018年5期)2018-11-06海上低渗储层防水锁强封堵钻井液技术钻井液与完井液(2018年5期)2018-02-13基于测井响应评价煤岩结构特征中国煤层气(2015年4期)2015-08-22随钻电阻率测井的固定探测深度合成方法西南石油大学学报(自然科学版)(2015年4期)2015-08-20高分辨率对地观测系统太空探索(2015年8期)2015-07-18基于Curvelet-Wavelet变换高分辨率遥感图像降噪浙江大学学报(工学版)(2015年1期)2015-03-01
    相关热词搜索: 测井 声波 成像

    • 文学百科
    • 故事大全
    • 优美句子
    • 范文
    • 美文
    • 散文
    • 小说文章