Using seismic surveys to investigate sediment distribution and to estimate burial fluxes of OC,N,and P in a canyon reservoir

Ningxiao Yu·Yong Qin·Feng Hao·Yunchao Lang·Fushun Wang

Abstract As a high-precision survey method, seismic surveying has been increasingly applied to inland water research,although its application to artificial reservoirs has remained limited. As a special artificial water body,reservoirs have important effects on the fluvial transport of material from land to ocean,and inevitably have complex terrain which can complicate and distort the results of seismic surveys.Therefore,there are still some problems need to be resolved in the application of seismic surveys in reservoirs with complex terrain.For this study,the Dongfeng Reservoir located in the upper reaches of the Wujiang River was chosen as an example to test the seismic survey method.Our testing showed that(1)because of the complex underwater terrain,the signal-to-noise ratio of the echo signal in canyon reservoir is low,making it difficult to determine sediment layers thicknesses in some areas;and(2)due to the large spatial heterogeneity of sediment distribution,insufficient density of cross-sections can lead to inaccurate interpolation results.To improve the accuracy of calculations,a mathematical method was used.Ultimately,the total burial mass of sediment was estimated at 2.85×107 tons, and the average burial rates of total organic carbon,total phosphorus,and total nitrogen were estimated at 0.194, 0.011, and 0.014 g cm-2 year-1,respectively.These values were close to the results of previous studies and hydrographic station data,indicating that seismic survey can be a reliable and efficient method for the mapping of reservoirs.

Keywords Dongfeng Reservoir·Seismic survey·

1 Introduction

As a link between land and sea,rivers play an important role in the global materials cycle(Sarmiento and Sundquist 1992;Vitousek et al.1997).It has been estimated that dissolved and particulate matter transported by rivers accounts for 90%of the flux from land to coastal seas(Liu et al.2009).For instance,about 2.7–2.9 Pg-C from land ecosystems enters inland water systems per year,where it experiences release,transformation,burial,and other processes during river transportation.Eventually,1.0 Pg-C(Regnier et al.2013)and approximately 48 Tg-N(Boyer et al.2006)enter coastal areas annually.The transport of nutrients(including C,N,P,Si,etc.)is significant for both oceanic and inland water ecosystems(Meybeck 1982),however,studies have shown that human activities,especially river damming,are affecting the role and function of rivers in the global materials cycle(Humborg et al.1997;Klaver et al.2007;Poff et al.2007;Regnier et al.2013).After river impounding, the slowed river flow and increased water surface area mean that hydraulic retention time is lengthened in reservoirs,and the fluidity of the water is limited.As a result,the‘‘river erosion effect’’,applicable in strong dynamic conditions, is gradually replaced by the‘‘lake sedimentation effect’’.Combined with the assimilation and absorption of aquatic organisms,some of the material carried by rivers can be trapped by dams and deposited in reservoirs(Hamilton and Schladow 1997;Salenc?on and The′bault 1996).Studies have shown that,globally,about 13%of the sediment transported by rivers becomes trapped in reservoirs annually(Beusen et al.2005).

Sediment deposited in reservoirs reduces their storage capacity and shortens the service life of dams(Bennett et al.2013);meanwhile,the deposition process leads to the accumulation in reservoirs of some substances(such as nutrients and heavy metals)which can cause potential environmental risks(Varol 2011).For instance,during the early diagenesis process,contaminants in sediments may be re-released in dissolved forms(e.g.dissolved phosphate),or in more toxic forms(e.g.methylmercury),into the overlying water(Balcom et al.2008;Pang et al.2008;Yin et al.2014).

Many studies have been conducted in order to estimate the material burial rates of reservoirs,but some significant unknowns remain.Traditional methods,such as radiometric dating(Crann et al.2015;Smoak and Swarzenski 2004),and sediment traps(Ritter 1972),can accurately measure the burial rate at discrete points.For most natural lakes with simple topography,hydrological conditions,and material sources,the burial rate and total buried amount of sediment can be roughly calculated for the whole lake by geometric models that combine the traditionally-collected data with parameters such as lake depth and area(Lehman 1975;Pajunen 2000).However,in reservoirs,due to the large spatial heterogeneity of sediment distribution,sediment accumulation rates can vary strongly between different sediment sampling sites(Mendonc?a et al.2014;Ferland et al.2012).For instance,in Lake Kariba,along the border between Zambia and Zimbabwe in Africa,sediment accumulation rates varied from 3.0×106to 3.0×107-t year-1across different areas(Kunz et al.2011).Therefore,it is inaccurate to estimate sediment burial rate for an entire reservoir based on a local area.

Recently, seismic sub-bottom profiling has become increasingly used in inland water research as a high-precision survey method(Mendonc?a et al.2014).Compared to traditional methods, seismic sub-bottom profiling can investigate the sediment distribution in detail.Due to the varying penetrability of acoustic p-waves in different media,sub-bottom stratum profiles with high resolution can be obtained by continuous cruising,allowing the total amount of sediment to be calculated(Wunderlich et al.2005).This technology has been studied and tested for decades(Zhu 2010)and is now starting to play an important role in geological exploration(Chen et al.2004),engineering guidance(Cui et al.2016),energy exploration,and other fields.However,seismic sub-bottom profiling is currently used less frequently for canyon-type reservoirs,and due to the complex underwater terrain of reservoirs it remains unclear(1)whether the underwater terrain of canyon-type reservoirs causes interference to acoustic signals;and(2)whether the results obtained by this method are reliable—and how the accuracy of results may be improved.

In this study,a dual-band(3.5 kHz and 10 kHz),seismic sub-bottom profiler was used to conduct a survey in the Dongfeng Reservoir(DFR),located in the upper reaches of the Wujiang River,in Guizhou Province,China.The main objectives were(1)to plot the isobaths and obtain sediment distribution maps for the DFR;(2)to evaluate and improve the accuracy of these maps;and(3)to calculate the accumulation rates of sediment,OC,N,and P in the DFR.

2 Study area and methods

2.1 Study area

The Wujiang River originates in Wumeng Mountain on the Yunnan–Guizhou Plateau,and is the largest tributary in the upper reaches of the Yangtze River.It has a total length of 1037 km and a drainage area of 88,267 km2(Wang et al.2010).The DFR is located in the upper reaches of the Wujiang River(106°05′–106°10′E,26°47′–26°51′N)and was built in 1994.It is a typical canyon reservoir,with seasonal regulation(Fig.1).The bedrock is mainly limestone and dolomite,and the river valley along this reservoir is typically flanked by cliffs.The DFR has a total storage capacity is 8.64×108m3, with a normal capacity of 4.91×108m3.The surface area under normal capacity is 19.06 km2,and the average hydraulic retention time is 33.4 days(Feng et al.2009).

2.2 Methods

2.2.1 Cruising survey and sampling

Low-frequency sound waves(＜1 kHz)have strong penetrability but low resolution of different geological strata,while high-frequency waves(＞10 kHz)have higher resolution but lower penetrability(Adams et al.2001;Hilbe et al.2011).For this reason,a dual-band(3.5 kHz and 10 kHz)seismic sub-bottom profiler(stratabox,SyQwest)was used in this study.The survey was conducted in August 2017.Cross-sections were taken with a separation distance of approximately 300–500 m,and point coordinates were recorded in real time by DGPS,then matched with sub-bottom profile data.During the seismic survey,sediment cores were collected using a gravity sediment core sampler in areas where the signals were good.

Fig.1 Map of sampling sites(HJD Hongjiadu Reservoir,PD Puding Reservoir,YZD Yinzidu Reservoir,red triangle hydrological station)

We chose two undisturbed sediment cores(from sample sites shown in Fig.1)for our chemical analyses.Because of the limitations of the sampling equipment,we were unable to collect complete sediment cores,and so the length of core D1 was 35 cm,while the length of core D2 was 43 cm.The cores were split at 1 cm intervals,and the divided samples were placed in clean,zip lock bags,which were then weighed.

2.2.2 Sub-bottom data processing and data fitting

Sub-bottom profile data were processed using Hypack software(Xylem,Inc.,USA)to obtain the longitude(X),latitude(Y),water depth(Z),and sediment thickness(H)at each observation point. Surfer (Golden Software, Inc.,USA)was used for the interpolation calculations to draw the isobath and sediment distribution maps.In calibrating our calculations,we assigned a value of zero to the sediment thickness and depth along the reservoir boundary.

In considering data interpolation,the Kriging method has better spatial correlation and can obtain more accurate interpolation results by editing the variogram r(h).However,the DFR layout is quite tortuous,the variograms for different regions varied significantly and in some cases were completely different,making it difficult to describe the whole reservoir with a single variogram.At the same time,the triangulated irregular network(TIN)method may generate unstable results for the interpolation of tributaries or strongly curved river reaches,and so was also not a suitable technique(Hua et al.2016).Due to these issues,the ‘‘inverse distance to power’’ method was used to interpolate the data in this study.

Considering the whole reservoir,the overall sedimentation mechanism is quite complex,which leads to a high spatial heterogeneity in the sediment distribution.In constrast,in smaller areas the hydrological conditions,water quality,and material sources could be regarded as sufficiently similar so that topography became the dominant factor affecting the sedimentation processes for these areas.On this basis,and according to the different hydrological conditions, the DFR was divided into several areas(Fig.1a–g),which were used to fit the relationship between the water depth and sediment thickness.According to the fitting results,we calculated the sediment thickness for areas where direct data could not be captured.

2.2.3 Chemical analysis

Sediment sample wet and dry weights were used to calculate sediment water content(WC).To measure total phosphorus(TP),samples were pretreated with perchloric acid and 98%sulfuric acid in a digestion tank for 12 h at 190°C,then measured by phosphorus molybdenum blue spectrophotometry (Jia et al. 2017). To measure total nitrogen(TN)and total organic carbon(TOC),samples were pretreated with 2 M hydrochloric acid to remove carbonates,and then the samples were washed to neutral with deionized water and freeze-dried.All samples were measured using an elemental analyzer.For particle size analysis,after the organic matter and the colloid were removed a dispersing agent(5%sodium hexametaphosphate solution)was added to the samples,and sizes were then determined using a Mastersizer 3000 laser particle size analyzer.

2.2.4 Calculation processes

In our calculations,we divided the DFR into many small grids,and the sediment thickness(H)of each grid in the sediment distribution map was obtained using Surfer.Because the physical and chemical parameters of the sediment were different at different burial depths,the sediment total mass was calculated as the sum of different sediment layers.Equations(1)and(2)were used for initial calculations,where Viis the volume of different stratification,? is sediment porosity,and sediment dry density ρsedis 2.5 g cm-3):

According to the elements’content at different sediment depths,the TOC mass(MTOC),TP mass(MTP),and TN mass(MTN)were calculated for each layer of sediment,using Eqs.(3)–(5):

Finally, we calculated the mass accumulation rate(MCR)of different elements using Eq.(6),in which 23 stands for the age of this reservoir,and S is its area:

3 Results and discussion

3.1 Variation of WC,TOC,TP and TN along the sediment profile

The content of clay(＜0.002 mm)and sand(0.02–2 mm)stayed relatively stable along the sediment profiles.The content of clay for sample D1 was between 30%and 40%,and between 25%and 30%for D2.The range of change of the sand content at both sites was about 10%.According to international soil classification standards,the sediments of the DFR fell within the categories of clayey silt generally,with some belonging to silty clay.Water content of surface sediment ranged up to 87.3%,decreasing gradually with depth, before stabilizing at approximately 50% below 30 cm(Fig.2).Due to equipment limitations,we were unable to measure water content at depths below 40 cm.

For deeper sediment,water content has been reported to be related to particle type and burial depth(Abraham et al.1999;Zhu et al.2017),so for the water content of sediment below 30 cm,we applied the results from Qin et al.(1983).

The TOC content in the surface sediment(0–5 cm)of reservoir site D1 was 3.05%,and 4.27%for the area in front of the dam(site D2),which was slightly higher than previously reported(Zhu 2005).Along the sediment profile,TOC showed a trend of decreasing concentration with depth,and stayed stable below 20–30 cm(Fig.3).TN showed a similar trend to TOC.In D1,TN content fell in a range between 0.14%and 0.29%;in D2,it ranged between 0.15%and 0.31%.TP variation was not as obvious as that of either TOC or TN. It ranged between 0.9 and 1.25 mg g-1in the two sites,which was basically consistent with previous reports(Yin et al.2010)(Fig.3).

The vertical distributions of these nutrients in the sediments was consistent with the diagenetic degradation threestage model for organic matter—namely,sedimentationdegradation-accumulation—described by Wang et al.(2000).The ratio of C and N is usually used to determine the source of organic matter in sediments.When the C:N ratio is ＜8,the organic carbon is mainly endogenous in origin;when the value is ＞12,the organic carbon is mainly exogenous(Ni et al.2011).The C:N ratios of surface sediments at both sites were all greater than 12,indicating that the recently deposited organic carbon is mainly of exogenous source.

Fig.2 Water content of the different sediment layers

In general,the content of nutrients in sediments can be more stable with increased deposit depth(He et al.2015).For instance,a study of the Ohio Reservoir has shown that the OC,N,and P content in sediments below 20 cm did not change noticeably (Vanni et al. 2011). Therefore, the average values of OC,N and P below 30 cm were used for calculations in this study(Table 1).

3.2 Characteristics of acoustic reflection signals in different representative areas

The lengths of sediment cores collected during the cruising survey were basically consistent with the sediment thicknesses obtained from the seismic signals,which indicated that the seismic sub-bottom profiler was reliable. The results showed an obvious spatial heterogeneity of sediment distribution.There was no sediment accumulation along the reservoir coast,but in deep and flat areas there were thick deposits(Fig.4a,b).This phenomenon indicated that the underwater terrain was an important factor in sediment deposition,as has been noted previously(Blais and Kalff 1995;Davis and Ford 1982).Two reasons have been suggested for this:(1)sediment in shallow water areas was easily resuspended during variations in water level;(2)as a canyon-type reservoir,the DFR has steep,near-shore slopes,and sediments easily slid into deeper areas under the combined influences of gravity and water flow(Blais and Kalff 1995;Davis and Ford 1982).

In addition,where the reservoir narrowed—especially at region B(Fig.1)—the sound signal was greatly scattered.The data signal-to-noise ratio,for both 3.5 kHz and 10 kHz bands,was low(Fig.4c),making it difficult to distinguish geological strata in these areas.Reasons for this probably include:(1)The complex terrain interfered with sound transmission;(2)The existing empirical formulas did not apply well to the complex underwater terrain,resulting in either over-compensation or under-compensation(Ji 2017).

3.3 Correction of the contour maps of water depth and sediment thickness

Before the reservoir’s construction,the depth of the original channel changed continuously.Figure 5a is the isobath map of the DFR drawn using the original data.However,the map indicated that the central channel water depth was discontinuous,which was inconsistent with the actual situation.In addition,the sediment distribution map lacked data in some areas(Fig.5b).Based on this map,a preliminary total sediments volume for the DFR was calculated as approximately 7.8×106m3,with 1.3×107tons total mass.These values were much smaller than the historical data from the hydrological station and previous studies,which indicated that at least 4.3×106tons of sediment was imported into the DFR annually over the period 1990–2003(Wu et al.2018).The main reasons for these differences may have been either defects of the interpolation method(Hua et al.2016),or missing data due to sound scattering. Therefore, we need a method to eliminate these errors and improve the accuracy of the results.

3.3.1 Isobath map of the DFR

There are many factors affecting material burial in a reservoir(Downing et al.2008).It is very difficult,and perhaps impossible,to describe the burial of reservoir material using a single variable.But,at a smaller scale,hydrological conditions,water quality,and material sources can all be taken to be similar,with topography the dominant factor affecting the sedimentation process for these areas,meaning that there could be some relationships between water depth and sediment thicknesses. Thus,obtaining more accurate water depth information for the DFR was the basis for resolving this relationship.

Fig.3 TOC,TN,and TP contents along the sediment profiles

Table 1 Physicochemical parameters of sediment cores

We imported the shape file of the original channel(before the reservoir was impounded)into the isobath map.Based on the measured values of different cross-sections,the original river channel was assigned a series of values(using the mean value of two adjacent sections),and then these values were added to the interpolation calculation to obtain a new isobath map(Fig.5c).

Fig.4 Characteristics of acoustic reflection signals in different representative areas

3.3.2 Data fitting and correction for the contour maps of sediment thickness

There was no significant correlation between water depth and sediment thickness in shallow areas.However,when the water depth was greater than a certain value(z′),sediment thickness was positively correlated with depth in a linear fashion(Fig.6),and this relationship was used in this study to fit the data for those regions where sediment thickness data was missing due to weak acoustic signals.We chose an appropriate z′(to maximize the r2of line L1),and values ＞z′were screened for linear fitting(Fig.6a).For areas with shallow water depth,line L1could not be used to calculate sediment thickness.Therefore,all data were fitted to a higher order polynomial to determine the curve L2.There was an intersection point(z′′)between the straight line L1and the curve L2(Fig.6).When the depth was＜z′′,curve L2was used for calculations,and when the depth was ＞z′′,L1was used for calculations.

Because of gravity and water flow,sediments in the DFR are mainly concentrated at the bottom of the central reservoir areas,and sedimentary deposits exhibited the characteristics of being small but deep.The straight line(L1)indicated the burial conditions of sediments after reaching a certain depth,so the correlation of L1was directly related to the accuracy of the fitting results.

Region A was close to the dam,with steep slopes on both sides—sediment thicknesses there showed good linear correlation with water depth(r2=0.92).Regions B and C belonged to the mainstream of the Yachi River.The relatively low correlation value(r2=0.57)of L1in region B could be caused by severe scattering of the sound waves,because the low signal-to-noise ratio made it difficult to determine sedimentary layers accurately.Region E was the confluence area of the Liuchong River and the Sancha River, and as such exhibited complex hydrology and topography.Consequently,the correlation of L1for this region was also relatively poor(r2=0.38).

The sediment input of the Sancha River is much less than that of Liuchong River,due to the construction of upstream reservoirs(Wu et al.2018).In region D(the Liuchong River),when the depth was ＞z′,linear correlation between sediment thickness and depth was also good(r2=0.64).Affected by the change in river width,the decreased slopes on both sides of regions G and H resulted in retention of some sediments on the underwater highlands,which led to relatively poor fitting results for the Sancha River.

Fig.5 Isobath map and sediment thickness map of the DFR(a isobath map before revise;b sediment distribution map before revise;c isobath map after revise;d sediments distribution map after revise)

3.3.3 Burial fluxes of sediment and nutrients

Using the revised water depth information,sediment thicknesses were calculated for local areas where data was not available.These fresh estimates were used to revise the sediment distribution map(Fig.5d),yielding new sediment volume and mass estimates of approximately 1.66×107m3and 2.85×107tons,respectively.

The Hongjiadu Hydrological Station is just downstream of the Hongjiadu Reservoir(HJD),which reflects the sediment imported from the Liuchong River into the DFR.But there are no long-term hydrological observation stations in the lower reaches of either the Puding or Yinzidu reservoirs.For this reason,sediment data from the Yangchang Hydrological Station,which is in the upper reaches of the Sancha River,was used in this study for comparative analysis(Fig.1).The annual sediment discharge recorded by the Hongjiadu Hydrological Station for the period 1990–2003 was approximately 4.33×107tons—which reduced to almost zero after 2004,due to the construction of the HJD Reservoir(Wu et al.2018).

The annual sediment flux gauged by the Yangchang Station from 1970–1983 was approximately 2.41×107tons.Then,due to the construction of the Puding and Yinzidu reservoirs in the lower reaches,large amounts of sediment were trapped in these reservoirs(Wu et al.2018;Xiang et al.2016).Therefore,sediments in the DFR were mainly from the Liuchong River during the period 1994–2004.The estimates calculated in this study were smaller than those based on the monitoring data of the hydrographic stations.Possible reasons for this discrepancy are(1)some sediment was deposited on the riverbed downstream of the hydrological station,and thus not transported all the way to the DFR;(2)after 1994,due to the construction of the upstream reservoirs(such as the HJD)and some soil–water conservation projects,sediment input to the DFR gradually reduced(Wu et al.2018);(3)impounding of the HJD Reservoir in 2004 led to a significant reduction in the amount of sediment transported by the Liuchong River.

In addition,the process of sand removal from the DFR has reduced sediment accumulation in some areas.Therefore,there was incomplete agreement between the DFR sediment accumulation estimates from this study and estimates based on data from the hydrological stations,although the results were basically consistent with the actual situation.

Fig.6 Fitting results for each region(red straight line:L1;black curve:L2)

Based on sediment accumulation data and chemical analyses,the burial fluxes of TOC,TN and TP in sediments were calculated (Table 2). These calculations yield a maximum MCROCfor the DFR of 2.317 g cm-2year-1,with an average value of 0.194 g cm-2year-1,which was much higher than the values reported for Hongfeng Lake(0.017 g cm-2year-1)and Baihua Lake(0.016 g cm-2-year-1),which are located in the Maotiao River(Wang 2001).However,these values also showed a positive correlation with the recharge coefficient(i.e.,the ratio of catchment area and lake area),which was consistent with previous studies(Downing et al.2008;Knoll et al.2014).The Wujiangdu Reservoir,located downstream and with a large recharge coefficient,has a much lower average mass accumulation rate of organic carbon(0.0071 g cm-2year-1,Yang et al.2017).This confirmed that upstream damming had obvious interception effects on material transport in the river,and that the relation between recharge coefficient and material burial rate was not applicable in cascade reservoirs.Our result was also consistent with previous estimates of organic carbon accumulation rates in temperate and cold zone reservoirs(up to 0.33 gC cm-2year-1,Mulholland and Elwood 1982),higher than the largest oligotrophic reservoir–Lake Kariba(0.0023 gC cm-2year-1,Kunz et al.2011),and lower than some small,eutrophic reservoirs(up to 1.74 gC cm-2year-1,Downing et al.2008).Values for MCRTPranged from 0 to 0.133 g cm-2year-1(with a mean value of 0.011 g cm-2year-1),and values for MCRTNranged from 0 to 0.183 g cm-2year-1(with a mean value of 0.014 g cm-2year-1).These results were similar to those reported for the Romania canyon reservoir,Iron Gate I(TP:0.015–0.033 g cm-2year-1, TN:0.07–0.355 g cm-2-year-1;Teodoru and Wehrli 2005).

Table 2 The sediment volume and MAR of different nutrients in the DFR

4 Conclusions

Seismic sub-bottom profiling is an efficient and low-cost survey method.By analyzing acoustic signals,geological information on the DFR was obtained,and the total volume of the sediments was calculated as approximately 1.66×107m3,with 2.85×107tons of total mass.Combined with sediment core collection and chemical analysis,the average mass accumulation rates of OC,TP,and TN in the DFR were estimated to be 0.194, 0.011 and 0.014 g cm-2year-1,respectively.

However,in canyon-type reservoirs,material burial is greatly affected by topography. Sediments are mainly distributed along the original channel of the reservoir area,exhibiting characteristically small areas with deep sediments.There are still some issues with seismic sub-bottom profiling,including(1)due to interference from terrain and other factors,the sound signal can become scattered in some areas;and(2)the sediment distribution area in canyon reservoirs is small,so that an insufficient density of cross-sections,or an uneven distribution of sample points can enlarge errors in the results of the data interpolation.

Because of the complex terrain in canyon-type reservoirs,increasing the transect densities and considering the original channel before planning the survey route can improve the accuracy of results. However, it is inevitable that sound scattering due to the terrain can still make it difficult to distinguish sediment layers in some areas,causing additional errors.

In this study,the sediment burial effect for the whole reservoir was evaluated using mathematical methods,with limited survey data.By collecting sediment samples,the contents in the sediment could be evaluated rapidly and accurately.Overall,we concluded that our study was a valuable contribution to future research on reservoir environments.

Acknowledgements We wish to thank our colleagues who contributed to this study:Huang Jiaqi,Jia Xiaobin,Zhou Tao,Wei Haobin and Chen Tianyu.They were instrumental in sampling and data processing.We acknowledge Li Shuang and Wai Wai Phyoe for their thoughtful comments on the manuscript.The study was funded by the National Natural Science Foundation of China(No.41573064)and the National Key Research and Development program of China(No.2016YFA0601003).