APP下载

近300年来福建长乐海岸沙丘记录的风沙环境演变*

2018-01-12于晓莉李志忠靳建辉赖海成申健玲苑秀全徐晓琳

海洋与湖沼 2017年2期
关键词:长乐风沙沙丘

于晓莉 李志忠① 靳建辉 赖海成 申健玲苑秀全 徐晓琳 程 延

(1. 福建师范大学地理科学学院 福州 350007; 2. 福建师范大学地理研究所 福州 350007;3. 湿润亚热带生态地理过程省部共建国家重点实验室 福州 350007)

海岸沙丘发生于大气、陆地和海洋交互作用的地带, 是记录海岸环境演变的良好信息载体(Murray-Wallace et al, 2001; Gardner et al, 2006; Scheffers et al,2012)。海岸沙丘发育与沙源供应、向岸盛行风、宽广的地形空间、趋于干冷的气候条件有关(王颖, 1987;吴正等, 1995a; 陈方等, 1996; 祁兴芬等, 2004; 刘建辉等, 2006)。近年来, 国内外学者对海岸沙丘的形成时代、发育过程、沉积构造及其记录的环境演变进行了广泛研究。Tamura等(2011, 2016) 通过探地雷达(GPR)探测和光释光(OSL)测年综合分析, 揭示了日本鸟取海岸横向沙脊在过去500年的演化过程, 并指出东亚冬季风控制了海岸风沙活动。Clemmensen等(2015)研究发现, 最近几百年来西北欧海岸大规模风沙活动主要发生在小冰期(LIA), 与北大西洋风暴增强呈正相关, 而同期 NAO值为连续的负值。Clarke等(2009)分析认为, 欧洲海岸和内陆(包括丹麦、苏格兰、威尔士等地)沙丘在小冰期风沙活动强烈, 是由于这一时期风暴潮活动增加, 促进了沙物质搬运和沉积。Dörschner等(2012)对台湾福隆海岸沙丘的研究表明, 海平面升降和构造运动是海岸沙丘不同发育阶段的主要控制因素, 大约 AD700前后的气候变化以及台风活动减少等促进了海岸沙丘的稳定。本文以我国亚热带海岸沙丘主要分布区福建闽江口以南的长乐海岸沙丘为研究对象, 通过系统的年代学分析并选用适宜的代用指标, 重建福建东南沿海小冰期以来的风沙活动序列, 为探索台湾海峡西岸小冰期海岸风沙活动特点及其成因寻找新的证据。

1 研究区概况

福建长乐为亚热带海洋性季风气候, 年平均气温19.3°C, 最热月均温27.6°C, 最冷月均温6—10°C。多年平均降水量 1382.3mm, 降水季节分配不均, 全年降水主要集中在 5—9月, 占全年降水量的 65%左右, 年相对湿度在 80%以上(谢皎如, 1998; 鹿世瑾,1999)。长乐地处台湾海峡西岸, 大风日数多, 夏季多台风, 全年≥8级的大风达30天左右。风向具有明显的季节性变化, 其中以 NE向占绝对优势, 成为研究区的主导风向, 仅 7月份以偏南风为主(图 1d), 年平均风速为 6.9m/s, 最大风速达 25m/s。长乐海岸线总体走向为NNE, 沿岸浪方向主要为SE和SW, 强浪方向为 ESE。沿海潮差大, 潮流强, 属于典型的强潮海岸(蔡爱智等, 2009)。长乐东部海岸风沙沉积物粒级组分主要是中砂和细砂, 海岸类型属于沙质海岸(陈方等, 1991)。

长乐东部毗邻福建最大河流闽江的入海河口南岸。全新世海侵后, 闽江流域泥沙逐渐充填了福州盆地及其峡谷口外浅水海湾, 钻探获知黏土和黏质粉砂层厚约2—6m。福州和长乐洼地的河口冲积层层厚5—10m, 主要是具平层理的泥和粉砂质泥。河口前积层的组成物是 2—4Φ的中细砂。河口浅滩区由两种地貌组成, 即河口水道向外延伸的槽沟和平行于槽沟不对称的河口坝。其中浅滩沙坝的不对称性主要受控于NE风引起的波浪, 表明了河口砂在优势风浪的作用下向南推移, 这一过程直接影响长乐海岸风沙地貌的发育和分布格局(石谦等, 1999; 刘苍字等,2001)。长乐海岸广泛分布典型的海岸沙丘, 在不同岸段的现代沙丘系统内侧, 风沙堆积形成了不同形态的早期沙丘系统, 包括横向沙丘链、穹形沙丘、抛物线沙丘、不规则沙丘等。其中长乐市江田镇东山村的沙丘类型为横向沙丘、横向沙丘链、横向复合型沙丘链(陈方等, 1992, 1996; 张文开等, 1995; 董玉祥,2006)。本研究所选取的沙丘类型是处在固定状态下的抛物线沙丘, 在 NE向盛行风作用下发育 NW-SE走向的沙丘脊。

2 材料与方法

2.1 沉积特征

本文采样工作于2015年7月8日在长乐市江田镇东山村东海沙丘剖面(DH)进行, 天气为阴天, 大风6—8级。DH剖面位于长乐市江田镇东山村, 省道201线南侧 50m 处(25°52′29.9″N, 119°34′44.5″E), 东部距离现代海岸线4.42km(图1)。东海村海岸沙丘以抛物线型沙丘为主, 开口大多朝向 NE, 两翼指向 NE或NNE方向, 大部分沙丘已遭到人为破坏。沙丘剖面在抛物线沙丘背风坡, 人工挖沙取土后形成, 整个剖面发育在一级海岸阶地上, 顶部已经完全被木麻黄林所覆盖, 出露厚度约6m, 底部海拔约19m, 未见底。从上到下剖面大致分为 3个层段(图 2a): A层: 0—130cm: 呈现灰白色(10YR7/2), 以中细砂为主, 其次为黏土和粉砂, 结构疏松, 缓倾斜层理, 倾向 NE、NNE, 顶部0—50cm含较多草本植物根系, 并含少量黑色植物碳屑。该层不整合的覆盖在下伏沙丘层上。B层: 130—500cm: 呈灰黄色(10YR7/4), 以中细砂为主, 其次为粗砂, 较紧密。其中 130—185cm 含较多植物根系和枯枝落叶, 发育缓倾斜层理, 倾角 2°—5°。185—500cm, 呈缓倾斜至水平层理, 其中灰白色砂层相间夹杂数条棕黄色砂层, 棕黄色砂层厚约2—5cm, 倾角约为 5°, 倾向 NW, 显示为埋藏的抛物线沙丘沉积构造。C层: 500—560cm, 呈淡棕褐色(2.5YR5/8), 发育波状层理, 较紧实, 湿度较大。

2.2 实验方法

2.2.1 OSL测年 根据剖面分层, 采用长度30cm、直径 6cm 的不锈钢管, 垂直打入去除表层浮土的剖面刻槽中采集测年样品, 取出样品后用黑色塑料袋和胶带密封钢管两端包装编号并带回实验室处理。共采集 6个 OSL样品, 从上往下依次编号为DH-01, DH-02, DH-03, DH-04, DH-05, DH-06。

图1 福建长乐东海海岸沙丘位置及区域风况Fig.1 Location and wind climate of DH (Donghai) coastal dune along the coast of Changle, Fujian

图2 福建长乐DH沙丘地层剖面及OSL样品采样位置Fig.2 Stratigraphic characteristics of the DH section in Changle, Fujian, and the optically stimulated luminescence (OSL) sampling location

在实验室前处理过程和样品制备如下: 在红光(波长为 640±10nm)条件下打开样品, 去掉钢管两端可能曝光的部分, 从这部分样品中取20g烘干测定环境剂量率和含水量, 并保留中心未污染、未曝光的样品进行等效剂量测定。将样品放入1000mL大烧杯中,依次加入 10%的盐酸(HCl)和 30%的双氧水(H2O2),除去样品中的碳酸盐和有机质。期间进行多次搅拌,待样品反应中止后加入氨水(NH3·H2O)进行中和, 然后用清水反复冲洗, 直至样品变为中性。之后将样品放入烘箱低温烘干(40°C), 筛选出<90μm的组分。将筛选出的组分用 30%的氟硅酸(H2SiF6)腐蚀以清除长石。再加10%的盐酸去除氟化物, 然后磁选去除磁性矿物, 用清水反复冲洗至中性后低温烘干, 最终筛选出58—65μm的石英颗粒, 制成样片, 待上机测验。

等效剂量(De)采用单片再生剂量法(SAR)进行测量, 流程如表1所示。样品的238U、232Th和40K含量由中国原子能科学研究院采用中子活化法(NAA)测得。根据Aitken(1998)提出的石英矿物的环境剂量率与238U,232Th,40K之间的转换关系, 并考虑各影响因素进行校正后得到环境剂量率(D)的值。OSL年龄计算公式为A=De/D, 单位为ka。

每个样品用 16个测片, 激发光源分别为蓝光二极管和红外二极管, 其中蓝光光源波长为 470±2nm, 半宽5nm, 最大功率为60mW/cm2; 红外光源波长为 830nm, 半宽 10nm, 测试中两种光源的最大功率为 90mW/cm2。本项测试采用两种光源的最大功率进行测量。光释光信号经由7.5mm厚的Hoya U-340滤光片进入 9235QA光电倍增管(PMT)进行记录, 人工辐照源为(90Sr/90Y)β源。所有样品使用福建师范大学Risø TL/OSL-DA-20C/D型全自动释光测量仪完成。

表1 SAR法测量流程(Murray等, 2000)Tab.1 Measurement process in SAR (single aliquot regenerative dose) method(Murray et al, 2000)

2.2.2 粒度分析 在DH剖面上, 按5cm等间距连续采集粒度分析样品91个。室内粒度测试前处理如下: (1)称取经过 2mm筛后混合均匀风干的待测样品5g, 放入洁净干燥的烧杯, 加入5mL10%的H2O2溶液静置 24h, 在电热板上加热充分反应以去除有机质;(2)加入5mL10%的HCl, 加热使其充分反应; (3)烧杯注满蒸馏水静置24h, 之后抽取上层清液, 重复3—5次直至溶液呈中性; (4)加入 5mL0.05mol/L的(NaPO3)6溶液进行分散, 并用超声波振荡3—5min完全分散之后上机测试。上机测试采用英国马尔文仪器公司生产的 MasterSizer2000激光粒度分析仪, 该仪器测量范围为 0.02—2000μm, 粒级分辨率为 0.01Φ,每个样品重复测试 3次, 以控制相对误差小于 2%,最后取其平均值。上述粒度实验在福建师范大学湿润亚热带生态地理过程省部共建国家重点实验室土壤粒度分析实验室完成。

3 结果分析

3.1 年代学

3.1.1 光释光测年分析 为了确保测试过程中测试条件的可靠性, 我们通过预热坪实验和剂量恢复实验, 确定了等效剂量的适宜预热温度为 220°C, 试验剂量的预热温度为180°C。沉积物样品OSL测年结果的可靠性很大程度上与样品的释光性质有关, 采用 SAR法一般要求样品的光释光信号以快组分为主(张家富等, 2009; Madsenet al, 2009)。从本文实验获得释光信号衰减曲线来看(图3), 样品天然OSL信号2s内快速衰减到本底值, 说明释光信号以快组分为主, 适合应用 SAR法进行等效剂量测定。再生剂量生长曲线采用线性拟合的方法, 可以看出各测量点基本都落在生长曲线上, 曲线拟合较好。从等效剂量值分布状况来看, 大多数样品等效剂量值比较集中,循环比率波动范围小于 5%, 表明测量过程中释光感量已经得到较好的校正。综上所述可以认为本文样品OSL测试数据基本可信。

从DH剖面6个样品的OSL测年结果(表2)可以看出, 剖面底部和顶部的年龄分别为302±9a和92±13a, 说明东海沙丘主要形成于近300a以来。并且各个层位样品的年龄值随深度的增加而变老, 基本符合正常层位的地层学关系(下老上新)。但是155cm处存在个别测年值倒置现象。从剖面分层特征看,155cm深度位于上下两段不同沉积特征的砂层中, 测片结果显示这个样品的年龄比较分散, 可能与样品的单测片测量结果分散、含水量估算误差和人为扰动作用等因素有关。

图3 福建长乐东海沙丘(DH)样品OSL信号衰减曲线、生长曲线和等效剂量分布Fig.3 Decay curves, dose response curves and equivalent dose in aliquots of quartz from DH section samples

表2 福建长乐东海剖面(DH)光释光测年结果Tab.2 The OSL results of DH section

3.1.2 独立年代学分析 为了验证 DH 海岸沙丘OSL测年结果的有效性, 需要对 OSL年龄进行独立年代学校验。通常运用海岸沙丘沉积中的贝壳碎片、古土壤层、泥炭进行14C定年(Aagaardet al, 2007;Madsenet al, 2009), 也可以根据环境考古方法校验沙丘发育的OSL年代。

根据前人研究, 闽江口海岸沙丘(福建长乐江田镇)下伏的褐黄色海相砂质黏土层14C测年为 2.47±0.09kaB.P.(程乾盛, 1989)。长乐、连江、宁德、霞浦等滨岸地带的海岸沙丘中富碳细砂层,14C年龄为1725±100aB.P(王雨灼, 1990)。福建南部漳浦的古雷半岛油沃村沙丘沙底部所含贝屑14C年龄为距今1400±90a(陈峰, 1986)。即福建沿海的海岸沙丘大部分是在2000年以来形成的。

前人在距离DH沙丘东北侧约300m处的东山沙丘(DS)剖面深度7m处采集的白色蚌壳AMS14C校正年代为 AD660左右, 在深度 7—7.4m大量出土的青釉瓷片鉴定为唐代中晚期制品(靳建辉等, 2015)。我们在 DH沙丘剖面 4—5m处也发现贝壳碎屑和陶瓷碎片, 虽然尚未进行测年实验, 但是根据其埋藏深度以及与DS沙丘的位置关系, 初步确定DH沙丘发育年龄晚于隋唐时期。DS沙丘的OSL测年显示其主要形成于 2kaB.P.以来(王贤立, 2015; 靳建辉等, 2015),而本文 DH沙丘的测年结果显示是在 302±9aB.P.(AD1713±9)以来发育的。可见, 这些毗邻海岸沙丘OSL测年界定的海岸风沙沉积时期, 与研究区含碳物质14C限定的起始风沙活动时间大体上是吻合的,总体上与华南海岸沙丘广泛发育时期相一致(吴正等,1995b), 同时可以看出历史上华南海岸沙丘的发育具有多期间断性特点。

3.2 粒度分析

粒度分析结果显示(图4), DH沙丘剖面的粒度组成以中砂和细砂为主, 两者总含量约 81.58%。其中,中砂含量 35.91%—60.74%, 平均含量 51.42%; 细砂含量12.38%—53.30%, 平均含量30.16%; 粗砂次之,平均含量为 14.96%; 极细砂含量较小, 平均含量为2.71%。黏土、粉砂和极粗砂含量很少, 三者总含量为 0.74%, 只在部分深度出现。例如黏土集中在 150—200cm部分, 粉砂在200cm以上部分显著增多, 极粗砂只出现在 200cm以上部分。我国海岸风成沙的平均粒径(Mz)介于 0.21—4.14Φ, 平均值为2.01Φ; 标准偏差(σ)介于 0.26—2.57Φ, 平均值为 0.6Φ; 偏度(Sk)介于–0.48—0.65, 平均值为 0.07; 峰态(Kg)介于0.80—3.09, 平均值为1.19。总体特征以细砂为主, 分选好,偏度对称近于正偏, 峰态很宽(董玉祥, 2003)。由图5可知, DH沙丘的平均粒径、标准偏差、偏度、峰态值分别介于1.30—2.20Φ, -0.97—0.56Φ, -1.82—1.03, 0.92—1.24, 平均值分别为1.71Φ, -0.66Φ, -0.19, 0.96; 总体特征以中细砂为主, 分选极好, 负偏至近对称和中等峰态。海岸风成沙的粒度组成一般以中砂或细砂为主,粉砂含量较少, 基本不含黏土(董玉祥, 2003; 王贤立,2015), 从整个剖面看, DH沙丘基本符合海岸风成沙的特征。

图4 福建长乐东海(DH)沙丘剖面不同粒级含量随深度的变化特征Fig.4 The variations of different grain size fractions with depths of the DH section

但是DH剖面200cm以上和460cm以下部分变化差异较大。在460cm以下部分, 极细砂和细砂含量显著增多。可能反映了上中下三个层段后期遭受的风化作用时间和强度不同。与北部约300m处的DS沙丘相比, DH沙丘的平均粒径偏粗, 可能与两个沙丘所处的海滩沙源、输沙通道位置有关, 但两者的平均粒径或中值粒径在剖面上的变化趋势大体相似。除了极粗砂, DH剖面样品粒度组分均有一定的波动旋回变化。因此, 本文采用DH沙丘平均粒径(Mz)作为记录沙丘发育时期风力强弱的替代指标。

图5 福建长乐东海(DH)沙丘粒度参数剖面变化Fig.5 Grain size characteristics of the DH section

从东海沙丘(DH)剖面代表样品的频率曲线(图 6)可以看出, 粒度频率分布区间大体一致, 粒度分布相对集中, 总体上呈现单峰宽度窄、峰值高的特点, 反映了相对稳定的沉积环境和较单一的物质来源。但是不整合面之上灰白色(10YR7/2)段与中部灰黄色(10YR7/4)段及下部棕褐色(2.5YR5/8)段粒度特征存在明显差异。上部黏土、极细砂、粉砂和极粗砂显著增多, 并呈现单峰加一细尾的特征, 细尾部分反映了海岸风沙沉积之后遭受了较强烈的化学风化作用,导致黏土和粉砂等含量增加。不整合面上下层样品的测年结果也说明了这一点。DH-01和 DH-02测年结果分别为 92±13a, 123±4a(图 2), 说明沙丘发育在123±4a(AD1892±4)以后可能经历了沉积间断, 之后再次活化接受沉积, 此期间可能经历了强烈的化学风化作用导致黏土和粉砂含量增加。中部层段除了不含有极粗砂, 其余各粒度组分变化不大。

下部棕褐色段极细砂、细砂显著增多, 中砂、粗砂和极粗砂减少, 频率曲线特征和中部表现一致。这种差异可能是海岸风沙沉积以后遭受后期化学风化作用造成的。岩性特征上表现为棕褐色, 类似闽东南沿海广泛分布的“老红砂”特征。吴正等(1995c)对华南沿海老红砂系统研究认为, 老红砂粒度组成以细砂和中砂为主, 含有粉砂和黏土, 这和 DH沙丘下部层段粒度特征一致。早期研究(吴正等, 1994; 曾从盛等, 1999)认为, 沙丘砂的红化作用时间要上千年至上万年, 才能形成淡红棕色。但近年研究发现, 在适宜气候条件下只要几百年时间, 海岸风沙就能发生红化现象。胡凡根等(2013)在福建晋江庵山沙丘研究中,就发现OSL年龄在1ka左右的老红砂, 本文剖面下部的“棕褐色”砂层年代很可能与晋江庵山老红砂沉积风化年代大体相当。

概率累积曲线能够反映沉积物的搬运方式并为解释沉积环境提供依据。从图7可以看出, DH沙丘剖面概率累积曲线分别呈现二段式、三段式和四段式的特点。0—130cm, 样品呈现滚动-跳跃-悬浮三段式特点, 其中滚动组分和悬浮组分直线斜率较大, 说明其分选性较好。130—285cm 层段, 样品呈现滚-跳-悬-悬四段式特点。粗截点在2Φ左右, 两个悬浮总体截点大约在8Φ。第一个悬浮总体区间宽, 斜率小, 分选差; 第二个悬浮总体区间窄, 斜率大, 分选性较好;285cm 以下, 基本呈现滚动-跳跃两段式特点, 不含悬浮组分, 其中跳跃组分百分含量较大, 达90%以上,滚动和跳跃组分直线斜率都比较大, 说明其分选性较好。一般说来, 风成沙丘沙具有跳跃组分高、滚动组分和悬浮组分含量低的特点(董玉祥等, 2008), DH沙丘剖面各层段跳跃组分含量都很高, 滚动组分和悬浮组分含量低, 因此是现代海岸典型的风沙沉积。另据萨胡判别函数:Y风成∶海滩= –3.5688Mz+3.7016σ2–2.0766Sk+3.1135Kg (Sahu, 1964)分析, 0—200cm样品Y风成∶海滩值介于–4.80—4.98, 属于风成沉积, 同时受到波浪、潮汐作用; 200—460cm样品Y风成∶海滩值介于–2.29—1.44, 以波浪、潮汐作用为主; 460cm以下样品Y风成∶海滩值介于–3.68—–1.07, 遭受风力、波浪和潮汐共同作用。

图6 长乐东海(DH)剖面代表样品频率曲线Fig.6 The grain size distribution frequency curves of the DH section

图7 长乐东海(DH)剖面各段样品概率累积曲线Fig.7 The probability accumulation curves of the DH section

4 讨论

4.1 长乐东海沙丘发育时期气候环境

关于华南海岸沙丘形成年龄的研究表明, 华南沿海现代海岸沙丘沙是近3000年以来的晚全新世发育的(吴正等, 1995b), 晚全新世是全球气候总体趋于降温且冷暖频繁交替的时期, 包括罗马暖期(约公元前1世纪—4世纪中期)、黑暗冷期(约4世纪末—10世纪前期)、中世纪暖期(约10世纪中期—13世纪末)、小冰期(15—19世纪)、20世纪暖期等(Lamb, 1977)。长乐DH沙丘发育于近300年来的AD1713—AD1923,即小冰期中晚期, 相当于晚清-民国时期。王绍武等(1998)利用史料、冰芯δ18O、树轮等重建的中国各区近 400—1000年的 10年平均温度序列表明, 1560—1699年和1790—1899年可以作为小冰期的两次主要寒冷阶段。

我国东部季风区近1000年温度序列和北部地区近 530年来的夏季气温序列表明, 小冰期内部存在1410—1499年、1560—1709年及1770—1919年3个冷期以及1500—1559年和1710—1769年2个相对偏暖的时期, 并且整体呈现冷干的气候特征(葛全胜等,2002; Liet al, 2007)。北半球温度变化曲线显示(图8),16世纪后期到17世纪(1570—1700年)、18世纪末至20世纪初(约1760—1920年)是过去500年北半球最明显的两个寒冷阶段, 18世纪气候则相对温暖。因此从年代学分析结果看, DH沙丘发育期是小冰期第二个冷期之后相对偏暖的时期和第三个冷期。

史料也记载了长乐东部海岸风沙地貌发育的气候背景。明弘治《长乐县志》记载: “长乐县东北有东湖, “延袤二十余里”, 乃由“海风飞沙, 积而成之”。清光绪二十八年(1902年)编《福建沿海图说》记载:“数十年前, 东海滨田舍相望, 今则一片平沙, 目断无人烟, 过其地者每不胜沧桑之感焉。”(林惠来,1982)。《平潭县志》记载: 平潭岛芦洋浦于乾隆十四年曾因大风, “海砂随潮涌上, 近海乡村悉遭压覆, 一夜砂埋芦洋浦十八村”。“横山东北, 清初该处共十三村, 雍正年间尽被风沙压废”。至1949年, 海岸林草殆尽, 一片荒芜。可以看出, DH沙丘发育阶段正是处于明清气候寒冷、海岸风沙活动强烈的时期。

4.2 小冰期长乐海岸沙丘发育的区域对比

如前所述, DH沙丘发育于近300年以来, 福建沿海其他区域也有小冰期(LIA)海岸沙丘发育的年代学记录。如长乐东山沙丘的 OSL测年结果记录了海岸沙丘发育的一个重要沉积期是在 800aB.P.左右(靳建辉等, 2015)。晋江庵山沙丘的 OSL测年结果显示沙丘是在距今854±95a形成的(胡凡根等, 2013)。漳浦霞美下蔡沙丘下覆泥炭层14C年龄为700±50aB.P.(陈承惠等, 1982)。这一时期福建气候总体表现为温凉背景下的干湿交替。例如闽北仙山盆地泥炭沉积物色度和有机碳含量变化指示了由偏湿(AD910—AD1640)至偏干(AD1640至今)的转化(宋瑞卿等, 2016)。

我国东部沿海其他地区也有小冰期海岸风沙沉积记录的研究报道。如辽宁盘锦海岸沙丘发育期为820±100 aB.P.以来, 发育初期冬季风较强, 且发育初期和末期气候均为干冷, 并对应北方降尘频率较高的时段(王小乐, 2014)。广东惠来南海乡下伏炭质粉细砂14C测年结果显示沙丘于510±90aB.P.开始发育(吴正等, 1995b)。台湾福隆海岸沙丘天然露头剖面的 OSL测年结果显示, 沙丘最近一次再活化是在AD1380—AD1840, 距该剖面1000m左右的四个沙脊测年结果显示沙丘发育在 AD1930—AD1970,AD1920—AD1960, AD1840—AD1940, AD1250—AD1560, 这一阶段沙丘稳定发育与过去 700年气候变化、沉积物供应增加、台风活动减少密切相关(Dörschner et al, 2012)。

全新世晚期尤其是 LIA时期, 世界海岸多个地区有海岸沙丘再活化和广泛发育的研究报道。例如东北亚的日本鸟取海岸(Tamura et al, 2011, 2016), 韩国江陵海岸(Choi et al, 2014), 南亚的印度那都海岸(Alappat et al, 2011), 欧洲西海岸和波罗地海岸等(Bailey et al, 2001; Clarke et al, 2002, 2006; Ballarini et al, 2003; Bateman et al, 2004; Wilson et al, 2004;Sommerville et al, 2007; Clemmensen et al, 2007,2015;), 这些地区的砂质海岸都曾处于小冰期强烈的冬季风作用下, 频繁的冬季风作用和风暴潮活动对海岸沙丘的发育产生重要影响。

4.3 长乐海岸沙丘发育与全球气候变化

小冰期成为世界海岸沙丘发育的一个活跃时期,与小冰期全球气候的剧烈波动有关(图8)。通过GPR探测、OSL测年和历史地图综合研究, Tamura等(2011,2016)重建了过去 500年日本鸟取海岸两条沙脊的演化过程, 发现东亚冬季风强度和海岸线变化对海岸沙丘发育有重要影响。例如, 15—17世纪和19世纪冬季风变强, 海岸沙丘向陆加积; 18世纪冬季风减弱,海岸沙丘向海进积。Van Vliet-Lanoë等(2014)将西北欧小冰期海岸沙丘发育与全新世风暴、北大西洋涛动(NAO)及大西洋年代际振荡(AMO)联系起来, 认为沙丘发育程度受控于北大西洋风暴强度。其中, 风暴增加的气候模式表现为NAO负相位、AMO正相位, 或NAO和AMO均处于负相位, 并在十年或更长尺度上表现出一定的相关性。

福建长乐海岸除了7月份为SSE风外, 全年其余月份盛行NE或NNE风。冬季NE向盛行风不仅风力大而且持续时间长, 直接控制了 DH沙丘的发育, 而这种NE风或NNE风受控于东亚冬季风系统。如果将DH沙丘记录的风沙活动期与广东湖光岩玛珥湖沉积物磁化率指示的冬季风强度(Yancheva et al, 2007)、格陵兰冰芯指示的北半球温度变化序列(Vinther et al,2006; 王绍武等, 2007)、北大西洋涛动指数NAO(Ortega et al, 2015)等对比(图8), 可以看出DH沙丘是在 LIA中后期亚洲冬季风强度稍弱时期发育的, 即AD1713—AD1923。此时, 北半球温度整体是偏低的。例如贵州荔波龙泉洞石笋记录显示, AD1800—AD1930我国西南地区的气候是寒冷-稍偏凉的(张美良等, 2009)。同时DH沙丘活跃期位于AD1650—AD1700和AD1820—AD1900两个低温期之一(图9)。因此DH沙丘发育期与广东玛珥湖沉积物磁化率指示的LIA以来东亚冬季风变化趋势、北半球温度变化序列有很好的对应关系。

Trouet等(2009, 2012)研究发现, LIA期间 NAO在正负相位之间反复振荡。当NAO持续为正值时, 冬季风暴增加, 促进西北欧海岸沙丘的发育; 当 NAO为负值时, 冬季反气旋活动会导致风暴减少, 抑制海岸沙丘发育。符淙斌等(2005)探讨了 NAO影响中国东部气候的大气动力学机制, 认为 NAO变化引起北大西洋环流系统的改变, 可能会导致西风带槽脊系统发生变化, 进而对地处欧亚大陆西风带下游的我国气候产生影响。在百年尺度上, 冬季NAOI(北大西洋涛动指数)及西伯利亚高压强度和东亚冬季风强度整体呈现明显的负相关关系, 相关系数分别达到–0.92、–0.95(王永波等, 2002)。DH沙丘发育于近300年左右, 这一时期气温是偏冷的, 并且 NAO与冬季气温的变化在AD1600年以后呈现明显的反相关关系(图9)。这一结论证实了NAO与冬季风强度在百年尺度上呈现反相关性这一观点。

福建地处我国东南沿海, 位于东亚夏季台风登陆的前沿区域, 海岸带受台风、风暴潮影响强烈。对过去1000年来长乐海岸沙丘(DS, DH)、晋江庵山沙丘的活动期与台风登陆福建沿海频率分析(Fogarty et al, 2006), 发现两者的相关性不明显。根据统计资料(1949—1990年)(杨华庭等, 1993), 福建沿海风暴潮最大增水 457cm发生在闽江口的梅花, 其余各站均在 200cm以下。由 DH沙丘距离现代东部海岸线约4.42km和DH剖面海拔高程12m可知, 现代风暴潮产生的风暴增水对其影响不大。但是分析DH剖面的沉积物粒度, 可推测明清时期风暴潮对其可能产生较大影响, 也可能是大风条件下, 强劲的 NE风吹蚀干燥的海滩砂堆积在 DH剖面。因此需要进一步提高海岸风沙沉积序列的OSL定年精度, 才有可能获得辨识台风、风暴潮事件较高分辨率的海岸风沙活动序列。

4.4 海岸沙丘发育与区域海岸线变迁

图8 福建长乐东海(DH)沙丘及其他海岸沙丘发育阶段与东亚冬季风强度、近千年北半球平均温度变化、北大西洋涛动指数(NAO)对比Fig.8 The timing of the DH coastal dune and other ones, and the comparison with the intensity of East Asian monsoon, Northern Hemisphere temperature, and North Atlantic Oscillation index (NAO)

海岸沙丘发育于沿海地带, 受海平面升降变化影响显著, 因而海岸沙丘的形成发育与海平面变化密切相关。海岸沙丘发育于高海面时期还是低海面时期, 目前仍存在不同认识, 主要集中在晚更新世古沙丘(老红砂)和全新世新沙丘。闽东南沿海的古沙丘(老红砂)主要是在晚更新世中期和晚期发育的, 当时大陆架广泛出露, 东海海面比现今低100m左右(曾从盛等, 1999), 全新世以来海平面则大幅上升。吴正等(1997)讨论了华南沿海老红砂的成因和形成环境, 认为老红砂是晚更新世末次冰期低海平面时形成的,并在其后高温湿热气候下和相对高海平面发生红化作用。胡凡根等(2012)对晋江科任老红砂的研究认为其发育于晚更新世末次间冰期以来气候暖湿的高海平面, 而在末次冰期低海面发生沉积间断。而全新世新沙丘(DS沙丘)记录了海岸风沙沉积大都发生在3000年以来福建海平面波动下降时期(王贤立,2015)。因此, 海岸沙丘发育时期的海平面高低与区域规模和时间尺度有关。

在晚全新世以来, 全球气候在温度波动下降中转为凉干, 世界各地普遍发生振荡性海退。福建闽江下游福州盆地 FZ5钻孔全新世孢粉和硅藻分析结果显示, 1900aB.P.以后福州盆地发生了快速海退(乐福远等, 2016), 至 700a B.P.以来海平面趋于稳定(谢在团等, 1983)。在中世纪(唐宋暖期)以后的 800aB.P.是长乐DS沙丘发育的一个重要沉积期, 可能记录了福建东部沿海小幅度海退的情况(王贤立, 2015)。距DH剖面北部约300m的DS剖面记录了海平面变化的信息, 其中0.8—0.2kaB.P.海平面是海退期、低海面时期,0.2kaB.P.至今, 海平面是海进期、接近现代海面(靳建辉等, 2015)。野外采样过程中, 在DH剖面淡棕褐色层段(5m)发现大量棕色瓷器碎片呈团块状分布, 可能是历史时期人类活动掩埋下来的。或许可以推断,DH沙丘发育时期海岸带存在一定的人类活动, 高海平面出现的可能性不大。

据历史文献记载(图 9), 福建长乐海岸线自秦汉以后, 经唐、宋、元、明、清、一直到现在, 总体上阶段性向海推进。这反映了晚全新世海平面总体下降背景下的次级波动情况(刘传标, 2014; 刘路, 2014,2015)。其中, 北线海澳(潭头澳)向海推进 1000—5000m 左右, 古海岸线之外的港口(福兴, 文石、霞江、汶上、凤洋等)现在都是潭头澳内的陆地。东北线, 现今的湖南镇, 过去的“沙洋坡”, 至今向海推进近万米。东南线, 鹤上洋, 漳港澳经围海造田, 至今向海岸推进约10000m。明弘治《长乐县志》记载:“长乐县东北有东湖, 至明初郑和下西洋时, 太平港虽开始淤塞变窄, 但主航道仍然深广, 所以郑和船队能在太平港内停泊。盖因这里地处闽江口南岸, 江水受东北风顶托, 把所带泥沙折向南流, 在沿岸形成沙堤和沙洲。到唐天宝年间, 鶠林调集乡民围筑而成湖泊; 相传, 这里的海岸线曾经三迁: 唐时岸线在今姚坑山、大当尾和龙口下, 清初移至海塘下, 清末又迁港嘴。”小冰期中后期AD1705—AD1920是长乐东部海岸沙丘的主要发育期, 至今长乐东海附近的海岸线向东推进约了2000—3000m。

图9 汉代以来各个历史时期长乐海岸线位置图(据文献《福建省历史地图集》修改)Fig.9 The coastline changes in Changle since the Han Dynasty in China (202—220BC)

5 结论

(1) 福建长乐海岸沙丘(DH)砂基本符合典型的海岸风成砂特征。砂样粒度组成以中砂和细砂为主,其次是粗砂, 三者总含量为96.54%。平均粒径、标准偏差、偏度、峰态值分别介于1.30—2.20Φ、–0.97—0.56Φ、–1.82—1.03、0.92—1.24, 平均值分别为1.71Φ、–0.66Φ、–0.19、0.96, 总体分选极好, 负偏至近对称和中等峰态。可能受风力作用强弱变化的影响,在剖面上各粒级组分和粒度参数均表现出一定的波动旋回变化, 因而可以作为过去风沙环境变化的有效代用指标。

(2) 根据OSL测年结果综合分析, 长乐东海沙丘(DH)主要是近 300年来的小冰期中后期 AD1713—AD1923发育的, 即晚清-民国早期是长乐海岸风沙沉积的活跃时期, 这一时期位于小冰期第二个冷期之后相对偏暖的时期和第三个冷期。虽然该时段东亚冬季风呈现波动减弱趋势, 但总体上呈现干冷多风的气候特点。大约AD1700以来, 在盛行东北风驱动下,长乐海岸风沙活动频繁, 叠加人类活动可能破坏了海岸带植被, 海岸带风沙地貌广泛发育, 长乐东海沙丘(DH)即形成于这一时期。

(3) LIA晚期以来, 福建长乐海岸沙丘发育期与东亚冬季风、北半球温度变化序列有很好的对应关系,即北半球温度偏低波动、东亚冬季风增强时期有利于海岸沙丘广泛发育。但是与北大西洋涛动(NAO)相位变化序列的相关性不明显。在较长的时间尺度内提取海岸风沙侵蚀和沉积记录, 可能要利用多个剖面加密采样进行 OSL定年确定年代学标尺, 以便更加准确地提取海岸沙丘记录的古气候信息。

(4) 长乐海岸沙丘(DH)发育与海平面变动密切相关。自汉代以来, 福建长乐海岸线逐渐后退, 由以前狭长的半岛加几个孤悬的海岛, 到唐代江海泥沙在闽江口淤积成浅滩低地, 宋代时继续向海推进发育为滨海平原, 明代发育为长乐两大平原直到广泛出现海岸风沙沉积, 海岸线向海推进了大约 5—10km。低海面及海退过程有效促进了研究区海岸沙丘向海进积, 同时也对明清以来人类经济活动产生较大影响。

王 颖, 朱大奎, 1987. 海岸沙丘成因的讨论. 中国沙漠, 7(3):29—40

王小乐, 2014. 盘锦海岸沙丘发育过程与环境指示意义. 北京: 中国地质大学(北京)硕士学位论文, 42—43

王永波, 施 能, 2002. 冬季北大西洋涛动与东亚大气环流及中国气候冷暖变化的多时间尺度相关. 见: 中国气象学会秘书处. 大气科学发展战略——中国气象学会第25次全国会员代表大会暨学术年会论文集. 北京: 气象出版社, 383—387

王雨灼, 1990. 福建省第四纪地层的划分. 福建地质, 9(4):289—306

王贤立, 2015. 福建长乐海岸沙丘记录的全新世海岸风沙环境. 福州: 福建师范大学硕士学位论文, 47—50

王绍武, 叶瑾琳, 龚道溢, 1998. 中国小冰期的气候. 第四纪研究,18(1): 54—64

王绍武, 闻新宇, 罗 勇等, 2007. 近千年中国温度序列的建立.科学通报, 52(8): 958—964

石 谦, 蔡爱智, 张金城, 1999. 闽江河口砂入海后的再搬运. 见:第九届全国海岸工程学术讨论会论文集. 北京: 海洋出版社,459—464

乐福远, 郑 卓, ROLETT B V等, 2016. 闽江下游全新世孢粉记录的植被与环境变化. 热带地理, 36(3): 417—426

刘 路, 2014. 海上丝绸之路起点福州甘棠港辨析(二). 福建社科情报, (6): 36—40

刘 路, 2015. 海上丝绸之路起点福州甘棠港辨析(三). 福建社科情报, (1): 35—40

刘传标, 2014. 海上丝绸之路起点福州甘棠港辨析. 福建社科情报,(5): 37—49

刘苍字, 贾海林, 陈祥锋, 2001. 闽江河口沉积结构与沉积作用.海洋与湖沼, 32(3): 177—184

刘建辉, 郭占荣, 2006. 福建长乐东部海岸沙丘发育成因及特征.福建地质, 25(4): 185—191

祁兴芬, 庄振业, 韩德亮等, 2004. 秦皇岛市海岸风成沙丘的研究.中国海洋大学学报, 34(4): 617—624

杨华庭, 田素珍, 叶 琳等, 1993. 中国海洋灾害四十年资料汇编(1949~1990). 北京: 海洋出版社, 98—170

吴 正, 王 为, 1997. 华南沿海老红砂的成因与古地理环境. 中国科学(D辑), 27(6): 537—542

吴 正, 黄 山, 金志敏等, 1994. 华南沿海老红砂的成因与红化作用. 地理学报, 49(4): 298—306

吴 正, 黄 山, 胡守真等, 1995a. 海岸沙丘的形成发育和形态类型. 见: 吴 正. 华南海岸风沙地貌研究. 北京: 科学出版社, 41—48

吴 正, 黄 山, 胡守真等, 1995b. 现代海岸沙丘的形成年龄.见: 吴 正. 华南海岸风沙地貌研究. 北京: 科学出版社,54—56

吴 正, 黄 山, 胡守真等, 1995c. 全新世环境变迁与新沙丘的发育. 见: 吴 正. 华南海岸风沙地貌研究. 北京: 科学出版社, 114—118

宋瑞卿, 朱 芸, 张江燕等, 2016. 福建仙山盆地泥炭沉积色度及其环境意义. 福建师范大学学报(自然科学版), 32(3):116—123

张文开, 李祖光, 1995. 福建长乐海岸沙丘形成发育及其区域分布特征. 中国沙漠, 15(1): 31—36

张美良, 朱晓燕, 程 海等, 2009. 贵州荔波 1200年来石笋高分辨率的古气候环境记录. 地球学报, 30(6): 831—840

张家富, 莫多闻, 夏正楷等, 2009. 沉积物的光释光测年和对沉积过程的指示意义. 第四纪研究, 29(1): 23—33

陈 方, 朱大奎, 1996. 闽江口海岸沙丘的形成与演化. 中国沙漠,16(3): 228—234

陈 方, 李祖光, 汪榕光等, 1992. 长乐东部沿海及海坛岛风沙地貌发育条件分析. 福建师范大学学报(自然科学版), 8(4):93—99

陈 方, 李祖光, 张文开等, 1991. 长乐东部沿岸风沙沉积物的粒度分布特征. 福建师范大学学报(自然科学版), 7(2): 84—91

陈 峰. 1986. 古雷半岛海滩岩的形成及闽南沿海海平面变化.见: 国际地质对比计划第 200号项目中国工作组编. 中国海平面变化. 北京: 海洋出版社, 166—172

陈承惠, 黄宝林, 王明亮, 1982. 闽南沿海全新世地质年代学研究.台湾海峡, 1(2): 64—73

林惠来, 1982. 台湾海峡西岸历史年代风沙的初探. 台湾海峡,1(2): 74—82

胡凡根, 李志忠, 靳建辉等, 2012. 福建晋江海岸带老红砂多期发育模式初步研究. 第四纪研究, 32(6): 1207—1220

胡凡根, 李志忠, 靳建辉等, 2013. 基于释光测年的福建晋江海岸沙丘粒度记录的风沙活动. 地理学报, 68(3): 343—356

符淙斌, 曾昭美, 2005. 最近530年冬季北大西洋涛动指数与中国东部夏季旱涝指数之联系. 科学通报, 50(14): 1512—1522

鹿世瑾, 1999. 福建气候. 北京: 气象出版社, 1—63

葛全胜, 郑景云, 方修琦等, 2002. 过去 2000年中国东部冬半年温度变化. 第四纪研究, 22(2): 166—173

董玉祥, 2003. 国内外海岸风成砂粒度参数特征的比较与分析. 中山大学学报(自然科学版), 42(4): 110—113

董玉祥, 2006. 中国海岸风沙地貌的类型及其分布规律. 海洋地质与第四纪地质, 26(4): 99—104

董玉祥, 马 骏, 黄德全, 2008. 福建长乐海岸横向前丘表面粒度分异研究. 沉积学报, 26(5): 813—819

程乾盛, 1989. 福建长乐全新世地层划分与对比的探讨. 福建地质,8(2): 100—109

曾从盛, 陈居成, 吴幼恭, 1999. 闽东南沿海老红砂沉积地层与形成时代. 中国沙漠, 19(4): 338—342

谢在团, 陈 峰, 刘维坤等, 1983. 福建全新世海滩岩与海平面变化. 台湾海峡, 2(1): 61—70

谢皎如, 1998. 闽东南海岸带季节性风沙气候环境. 福建地理,13(1): 6—11

靳建辉, 李志忠, 胡凡根等, 2015. 全新世中晚期福建海岸沙丘记录的海岸环境与人类活动. 地理学报, 70(5): 751—765

蔡爱智, 石 谦, 2009. 台湾海峡成因初探. 厦门: 厦门大学出版社, 26—32

Aagaard T, Orford J, Murray A S, 2007. Environmental controls on coastal dune formation; Skallingen Spit, Denmark.Geomorphology, 83(1—2): 29—47

Aitken M J, 1998. An Introduction to Optical Dating: The Dating of Quaternary Sediments by the Use of Photon-Stimulated Luminescence. London: Oxford University Press

Alappat L, Frechen M, Ramesh R et al, 2011. Evolution of late

Holocene coastal dunes in the Cauvery delta region of Tamil

Nadu, India. Journal of Asian Earth Sciences, 42(3): 381—397

Bailey S D, Wintle A G, Duller G A T et al, 2001. Sand deposition during the last millennium at Aberffraw, Anglesey, North Wales as determined by OSL dating of quartz. Quaternary Science Reviews, 20(5—9): 701—704

Ballarini M, Wallinga J, Murray A S et al, 2003. Optical dating of young coastal dunes on a decadal time scale. Quaternary Science Reviews, 22(10—13), 1011—1017

Bateman M D, Godby S P, 2004. Late-Holocene inland dune activity in the UK: a case study from Breckland, East Anglia. The Holocene, 14(4): 579—588

Choi K H, Choi J H, Kim J W, 2014. Reconstruction of Holocene coastal progradation on the east coast of Korea based on OSL dating and GPR surveys of beach-foredune ridges. The Holocene, 24(1): 24—34

Clarke M, Rendell H, Tastet J P et al, 2002. Late-Holocene sand invasion and North Atlantic storminess along the Aquitaine coast, southwest France. The Holocene, 12(2): 231—238

Clarke M L, Rendell H M, 2006. Effects of storminess, sand supply and the North Atlantic oscillation on sand invasion and coastal dune accretion in western Portugal. The Holocene, 16(3):341—355

Clarke M L, Rendell H M, 2009. The impact of North Atlantic storminess on western European coasts: a review. Quaternary International, 195(1—2): 31—41

Clemmensen L B, Bjørnsen M, Murray A et al, 2007. Formation of Aeolian dunes on Anholt, Denmark since AD 1560: a record of deforestation and increased storminess. Sedimentary Geology,199(3—4): 171—187

Clemmensen L B, Glad A C, Hansen K W T et al, 2015. Episodes of Aeolian sand movement on a large spit system (Skagen Odde,Denmark) and North Atlantic storminess during the Little Ice Age. Bulletin of the Geological Society of Denmark, 63: 17—28

Dörschner N, Reimann T, Wenske D et al, 2012. Reconstruction of the Holocene coastal development at Fulong Beach in north-eastern Taiwan using optically stimulated luminescence(OSL) dating. Quaternary International, 263: 3—13

Fogarty E A, Elsner J B, Jagger T H et al, 2006. Variations in typhoon landfalls over China. Advances in Aomospheric Sciences, 23(5): 665—677

Gardner T W, Webb J, Davis A G et al, 2006. Late Pleistocene landscape response to climate change: eolian and alluvial fan deposition, Cape Liptrap, southeastern Australia. Quaternary Science Reviews, 25(13—14): 1552—1569

Lamb H H, 1977. Climate: Present, Past and Future. Volume 2:Climatic History and the Future. London, England: Methuen and Co. Ltd., 835

Li J B, Chen F H, Cook E R et al, 2007. Drought reconstruction for North Central China from tree rings: the value of the Palmer drought severity index. International Journal of Climatology,27(7): 903—909

Madsen A T, Murray A S, 2009. Optically stimulated luminescence dating of young sediments: a review. Geomorphology,109(1—2): 3—16

Murray A S, Wintle A G, 2000. Luminescence dating of quartz using an improved single-aliquot regenerative-dose protocol.Radiation Measurements, 32(1): 57—73

Murray-Wallace C V, Brooke B P, Cann J H et al, 2001. Whole-rock aminostratigraphy of the Coorong Coastal Plain, South Australia:towards a 1 million year record of sea-level highstands. Journal of the Geological Society, 158(1): 111—124

Ortega P, Lehner F, Swingedouw D et al, 2015. A model-tested North Atlantic Oscillation reconstruction for the past millennium.Nature, 523(7558): 71—74

Sahu B K, 1964. Depositional mechanisms from the size analysis of clastic sediments. Journal of Sedimentary Research, 34(1):73—83

Scheffers A, Engel M, Scheffers S et al, 2012. Beach ridge systems-archives for Holocene coastal events?. Progress in Physical Geography, 36(1): 5—37

Sommerville A A, Hansom J D, Housley R A et al, 2007. Optically stimulated luminescence (OSL) dating of coastal aeolian sand accumulation in Sanday, Orkney Islands, Scotland. The Holocene, 17(5): 627—637

Tamura T, Bateman M D, Kodama Y et al, 2011. Building of shore-oblique transverse dune ridges revealed by groundpenetrating radar and optical dating over the last 500 years on Tottori coast, Japan Sea. Geomorphology, 132(3—4): 153—166

Tamura T, Kodama Y, Bateman M D et al, 2016. Late Holocene Aeolian sedimentation in the Tottori coastal dune field, Japan Sea, affected by the East Asian Winter monsoon. Quaternary International, 397: 147—158

Trouet V, Esper J, Graham N E et al, 2009. Persistent positive North Atlantic Oscillation mode dominated the Medieval Climate Anomaly. Science, 324(5923): 78—80

Trouet V, Scourse J D, Raible C C, 2012. North Atlantic storminess and Atlantic Meridional Overturning Circulation during the last Millennium: reconciling contradictory proxy records of NAO variability. Global and Planetary Change, 84—85: 48—55

Van Vliet-Lanoë Van B, Penaud A, Hénaff A et al, 2014. Middle- to late-Holocene storminess in Brittany (NW France): Part II-the chronology of events and climate forcing. The Holocene, 24(4):434—453

Vinther B M, Clausen H B, Johnsen S J et al, 2006. A synchronized dating of three Greenland ice cores throughout the Holocene.Journal of Geophysical Research: Atmospheres, 111(D13):D13102

Wilson P, McGourty J, Bateman M D, 2004. Mid-to late-Holocene coastal dune stratigraphy for the north coast of Northern Ireland.The Holocene, 14(3): 406—416

Yancheva G, Nowaczyk N R, Mingram J et al, 2007. Influence of the intertropical convergence zone on the East Asian monsoon.Nature, 445(7123): 74—77

猜你喜欢

长乐风沙沙丘
沙丘
长乐姑娘陈佼怡 新晋欧洲华姐冠军
一片冰心育桃李 六秩弦歌谱华章——福建省长乐华侨中学简介
会弹琴的沙丘
时间的年轮
东明县风沙化土地监测与治理
都怪祖先
沙丘
大宋小神探 夜黑黑,敲窗声
国庆记忆