AN IMPROVED DEM CONSTRUCTION METHOD FOR MUDFLATS BASED ON BJ-1 SMALL SATELLITE IMAGES : A CASE STUDY ON BOHAI BAY

The topographic measurement of muddy tidal flat is restricted by the difficulty of access to the complex, wide-range and dynamic tidal conditions. Then the waterline detection method (WDM) has the potential to investigate the morph-dynamics quantitatively by utilizing large archives of satellite images. The study explores the potential for using WDM with BJ-1 small satellite images to construct a digital elevation model (DEM) of a wide and grading mudflat. Three major conclusions of the study are as follows: (1) A new intelligent correlating model of waterline detection considering different tidal stages and local geographic conditions was explored. With this correlative algorithm waterline detection model, a series of waterlines were extracted from multi-temporal remotely sensing images collected over the period of a year. The model proved to detect waterlines more efficiently and exactly. (2)The spatial structure of elevation superimposing on the points of waterlines was firstly constructed and a more accurate hydrodynamic ocean tide grid model was used. By the newly constructed abnormal hydrology evaluation model, a more reasonable and reliable set of waterline points was acquired to construct a smoother TIN and GRID DEM. (3) DEM maps of Bohai Bay, with a spatial resolution of about 30 m and height accuracy of about 0.35m considering LiDAR and 0.19m considering RTK surveying were constructed over an area of about 266km. Results show that remote sensing research in extremely turbid estuaries and tidal areas is possible and is an effective tool for monitoring the tidal flats.


INTRODUCTION
The waterline detection method (WDM) based on satellite images is one of the most effective methods for constructing digital elevation model (DEM) for tidal flats (Liu Y.X. 2012).Typical research included Mason et al. (Mason, 2001, 2010).Waterlines were extracted from satellite Synthetic Aperture Radar (SAR) images in the tidal flats of Morecambe Bay from 1991-1994.Waterlines were heighted using a hydrodynamic model and tide gauge data and then constructed the intertidal DEM.Won, J.S. et al. applied ERS-1/2 tandem pairs and spaceborne radar interferometry (InSAR) to the Korean tidal flats to acquire waterlines, then the tidal flat DEM were constructed and a vertical accuracy of higher than 20 cm was required (Won, J.S.et al. 2003a(Won, J.S.et al. , 2003b)).In the study of Damien et al.,(2009) waterline method was used to build a high-resolution DEM on a beach of Aber Benoit(France) by recording the shoreline movements during a whole tide rise with a thermal infrared camera, and with simultaneous check of the water height(DGPS on board a boat) to provide an absolute calibration.The above methods treat the waterline as an altimeter that change with water level in the intertidal zones, but factors including tidal delay effect and coastal terrain are not considered.
Chinese researchers have applied multi-temporal remotely sensed images such as Landsat TM, ETM and SPOT to retrieve DEM of the muddy tidal flat.Wang Y.J. (2003) built the correlation between soil water content and topography of intertidal flat through statistics by Landsat TM images and simulated the topography of intertidal flat in DaFeng, Jiangsu province.Zheng Z.S.et al. (2007) used waterline detection methods and hydrodynamic model Delft3D to assign elevation to the waterlines at the satellite overpass in Dongtan of the Yangtze estuary.Shen F. et al. (2008) performed an approach to construct DEM by a series of waterlines extracted from multitemporal Landsat TM data.The elevation of each waterline was assigned into tidal level value by the estimation from tidal station record.He M.B. et al. (2008) aimed at the wide Jiuduansha intertidal mudflat and used water points instead of waterline to acquire elevation with RTK-GPS.These researchers in China have all used images aboard the country, so it will be restricted in acquiring images with proper time and good quality.And furthermore, the waterline elevation models haven't applied the high resolution tidal grid model of China near-shoal coast.
In this study, firstly and adequately considered the tidal current transmitting pattern of China near-shoal, a new intelligent correlating model of waterline detection considering different tidal stages and local geographic conditions was explored.By using high-resolution 1.2′×1.2′Chinanear-shoal tide model, the spatial structure of elevation superimposing on the points of waterlines was firstly constructed and a more reasonable and reliable set of waterline points was acquired to construct a smoother TIN and GRID DEM.

Estuarine Delta
Bohai Bay is one of the large bays of China, and it belongs to the typical mudflat plain coast.The study area located in this bay with the coordinates: 38°20′1″～38°42′12″North Latitude and 117°31′32″ ～ 117°52′22″ East Longitude.We conducted the experiment on Bohai Bay, from Tanggu Harbor to Huanghua Harbor.The intertidal zone along Bohai Bay coast is very wide and silty, where the rivers carry magnificent mud and sand.The main rivers are Yellow River, Haihe River, Jiyunhe Canal and Luanhe River.The experiment site is among the Luanhe River estuarine delta, Yongding River water system plain and Yellow River estuarine delta.The area extends from the edge of the plain, which is mainly made up of sand and silt and has tiny slope.Straight coastline, tiny slope and wide mudflat are the characters of the coastline.The area develops from dominant Yellow River water system and is at the edge of the Huabei Plain.The main slope is very small and the inlet material is tiny sand and clay.The area forms typical tidal mudflat plain coastal zone, which is dominated both by the dynamic tidal effect and the deposit of silt and sand (Figure1).

Scope and Area
In 1980, the mudflats coast spread the length for about 641.7km in the whole Bohai Bay from Liaohe plain to Huabei plain.But in 2010, mudflats coast was only 158km left, and the mudflat coast with the length of 483.7km had been turned into breed aquatics or fishery land or harbor construction area.The study area was a promising region for exploitation of marine resources and harbor construction.Between the year 2009 and 2012, the amount of mudflats were exploited.Typical harbor constructions include engineering of Tanggu Harbor, Dagang Harbor in the south of Tanggu and Qikou in Hebei Province and Huanghua Coal Harbor with first period.So the mudflat DEM construction is with great importance to the coastal zone exploitation and sea route silt cleaning.

Slope and Zonal area of tidal flat
The width of tidal flat near Huanghua Harbor is about 15-16 kilometers.Therein the width of intertidal flat is about 5-6 kilometers, and the width of lower tidal area is about 9-10 kilometers.The mean slope of the intertidal flat is about 0.6‰ and the slope of lower tidal area is about 0.44‰.The tidal flat develops in an open sea area, and it has many differences from the closed sea area.The substance of the intertidal flat is mainly mud and muddy sand (Wang F, 2011).The mean tidal range varies from 2m to 4m, and the local maximum tidal range (Tanggu Gauge) can reach 5.1m (Li X.B., 2011).

Maps and charts of this area
For the intertidal zones have been mapping blank areas for years, terrain maps, coastline and depth ordnance survey maps and Hydrographic Office charts of Bohai Bay intertidal zone are all very old, last updated 1983(some area last updated 1959).The chart of Qikou and the adjacent Area (No.5208) with the scale 1:50 000 produced in 2008, but the adopted data were still from the old maps.

MEASUREMENTS AND DATA
The range of tidal flat is very wide and with very tiny slope.Waterlines in such area will be well distinguished with 30 meters' spatial resolution.Considering the economy and history sequences, three kinds of middle spatial resolution multispectral images were selected such as TM/ETM+、HJ-1 and BJ-1 small satellite images.

Analysis of acquiring data
Based on a series of China property Beijing-1(BJ-1) small satellite image data, which were kept in achieves and can cover an area of 600km*600km, a series of images were selected with all the tidal conditions of the mud tidal land.In a year of spring, summer, autumn and winter, images with sunny weather were selected with the temporal sequences, such as 2009-03-06, 2009-04-26, 2009-06-29, and 2009-10-17.The tidal conditions included high tide, middle tide and low tide.

Extraction of waterline
The delineation of waterlines of optical remote sensing data of tidal flat areas is hampered by optical complexity and often extreme turbidity.As the land-sea boundary has always be mixed up with suspended mud and sand, the remnants water on the flat and the wind and wave, the fuzzy degree of image and the high, middle or low tide situation of the waterlines will be very different.The correlation model of waterline detection and the accuracy assessment models were put forward aiming at the four kinds of waterlines (Table .1).Each kind of waterlines could be efficiently and exactly detected and assessed.An object-based segmentation method was applied in detecting fuzzy or clear waterline of high tide with improved accuracy.Details are given in (Wu D. et al., 2014).As for fuzzy images of low tide, a newly created method considering BJ-1 small satellite integral time called water index was applied.Different form NDWI specially for TM images Band 2(Green)and Band 4(MIR), the new water index considered BJ-1 images Band 1(Green)and Band 3(NIR), and a threshold was added correlating BJ-1 integral time which was correlated with gray scale of muddy suspended water body (Wu D. et al., 2014).And for clear images of low tide, Canny arithmetic operators had good detection accuracy as edge an improved detection method.

Table.1 Analysis of BJ-1 images and Waterline detection method
A series of waterlines were extracted from 30 multi-temporal remotely sensing images in Bohai Bay over the period of a year using the above approach.This correlation model rose up the efficiency and reduced the complexity of the tidal flat conditions.Further the waterlines should be with high assurance.

Improved waterline heighting
Waterlines detected from remotely sensed images reflect the instantaneous land and sea interface while the satellite is passing through.The waterline is commonly considered as a contour.Then the geocoded positions of the waterline are superimposed the heights relative to the tide gauge data nearby.
If it is lack of in situ data, estimation data from tidal gauge station record could substitute.But according to the theory of the marine tide, the promulgating of tide has certain character and discipline, so the points on the waterline would not be with a same elevation, and the whole waterline could not be considered as a contour.In an attempt to achieve the vertical accuracy of the waterline heighting, three approaches were studied.Firstly we built the waterline elevation spatial structure, then the fined tidal grid model was adopted, and thirdly, the waterline point's abnormal model was set up to improve the waterline heighting accuracy.

Elevation spatial structure
Waterline records the up and down of the instantaneous sea level, then the height of the points on the waterline is relative to the local tide character.In marine sounding, the tide gauge station will be set to gain the local tide character.With a long period of tide gauge, the local tidal assimilation constant and the instantaneous water level correction value of different temporal will be achieved.The water level correction idea was applied into the calculation of waterline elevation in this paper.
A spatial structure of waterline (water point) elevation was constructed based on the marine sounding spatial structure.
Through the conversion of sounding and elevation, the water point elevation will be computed.

Figure.2 Spatial structure of waterline elevation
From Figure .2, the relation will be built up between the instant water level field and the instant waterline elevation field; the elevation h will be calculated by the following equation, is not related with time, so it has stability.As for spatial conversion, h and  are all part of T ， h can be calculated by  and T .Thus, the substantial conversion is the spatial structure of the sounding datum and the elevation datum.The key is the calculation of the instant water level ) , , ( t y x T .

High resolution tidal grid
To acquire ) , , ( t y x T , the tidal characters should be studied in the research area.Shown as the waterline elevation spatial structure model, the most reliable acquiring method will be the in situ tidal gauge.But the actual work may not allow setting tide gauge station in any near shoal position or in any density.The limited and discrete tide gauge stations are used to build a zonal interpolated model, supposing the tidal level along the coast is evenly promulgated.And this model is fit for the tidal area where the Co-pulse lag line of main tidal component is vertical to the coast.However, the instantaneous sea level could not be simply calculated by the zonal interpolated model.As showed in Figure.Based on the POM model and "blending" assimilation, the altimetry data analysis results are assimilated, and the finer 1.2′×1.2′regional tidal grid model of the China sallow water areas is built up, the standard deviation of the 9 main tidal components is 12.5cm(Zhang J. 2011).
For a small scope of area or an inlet of a small river, the waterline can be directly heighted with water level of the tide gauge station.But for a middle scale area where the Co-pulse lag lines are vertical to the coastline, the zonal interpolated model will be used to correct the direct evaluation error.For the large scale area where the coastal shallow sea water has sophisticated tide such as the amplitude and the pulse lag are not in-phase, the high resolution regional tide grid model will be adopted.We implemented the three models respectively in Bohai Bay coastal shallow water area, the standard deviations of the directly heighted method, the zonal interpolated model and the 1.2′×1.2′tide grid model were 30cm, 20cm and 12.5cm.

Waterline points abnormal
Study the waterlines set, contradicted phenomena appeared to the distribution of the waterlines.Sometimes they were mixed together.One point might appear simultaneously on different waterlines (Figure .5).We call this phenomenon the waterline points abnormal.It was caused by various hydrologies abnormal, such as sea condition, sea wind, storm tide, cyclone etc.

Figure.5 Waterlines set on the tidal flat of BJ-1 images
In order to correct the deviation of the waterline points abnormal, a qualitative and quantitative evaluation model was put forward.
We conducted an in-situ transect surveying, combined with the geography and meteorology background in the study area.Index numbers of each waterline were ranked from high level to low level at the tidal flat.Then a weight assignment was imposed on each waterline with formula (2).  Figure .7 the optimal discrete waterline points

Construction of DEM
The set of points with accurate elevation values acquired with the above model is used to construct Triangulated Irregular Networks (TINs).Then an interpolation for each grid elevation was performed in accordance with the associated triangle.
In these cases of Bohai Bay, north east of China, shoreline information extracted from multi-spectral images of BJ-1 small satellite at different tidal levels and the optimal waterline points were interpolated into a digital elevation model.DEM maps with a spatial resolution of 30 meters were constructed over an area of about 266km 2 .The North aspect is 0°and calculating in clockwise, the value scope is 0°~360°.In ArcGIS software, there are 9 aspects: tiny slope(-1), north slope(0°2 2.5°，337.5°~360°),north-east slope(22.5°~67.5°),east slope(67.5°~112.5 °), south-east slope (112.5 °157.5 °), south slope(157.5°~202.5 °), south-west slope(202.5°~247.5 °), west slope(247.5°~292.5 °) and north-west slope(292.5°~337.5°).As for the results, lots of the aspects of the gradient were eastern and north-eastern.So the overall slope leant from land to sea and from south to north.The geographical explain was that the phenomena resulted in the tidal dynamical character and tidal flat aggradations character.Furthermore the aspects were caused by Qikou hydrology transect and Bohai Bay morphology sediment and erosion.The mechanism of the geography showed that the mudflat was mainly influenced by the sand transportation from the southern Huanghe River (Wang Y., 2012).So the aspect leant from south to north.This study not only gave a quantitative prove to the geographic discipline, but also filled up the mapping gap for this area.

LiDAR validation
For validation, LiDAR data from National Ocean Office were used which acquired on 24, June, 2009.The resolution of the result grid was 2 meter ×2 meter.For the inversion DEM grid resolution was 30 meter ×30 meter, we resampled the LiDAR data to the same resolution.The LiDAR DEM chart (30m×30m) and the subtract error distribution map of the inversion DEM and LiDAR DEM were showed as Figure.

In-situ RTK validation
We also performed ground leveling in order to validate the inversion precision.From the latest monitoring BJ-1 image, three profiles were selected that had no obvious manmade construction and inning since 2009.The profile was transected the gently sloping tidal flat spreading from coastline to semitidal, the total surveying length was 5.2 kilometers.The transects adopted RTK detailed surveying, with the horizontal positioning precision 10mm+1.0×10 - D and vertical positioning precision 20mm+1.0×10 - D (D refers to the distance between the reference station and the mobile station ).There were 178 detailed points, and the interval was not definitely equal distributes, the intervals varied from 20 meters to 100 meters.Total ground leveling results were all RTK fixed solution.For further quantitative evaluation of the inversion precision, the absolute error was calculated by inversion elevation subtracts surveying elevation, and the relative error was calculated by inversion elevation divides surveying elevation.
The average precision of the waterline points was 100% subtracts the average relative error.As showed in Table .2,statistic errors were given including the maximum, minimum, average and RMS (Root Mean Square) of the data set in the three transects.The elevation average absolute error was 0.1880m, or 18.8cm (19cm), the average relative error is 11.3% (Table .2).The average precision of elevations of the waterline points set is 88.7%.From the transect validation result, the waterline method in the paper has reached a good inversion precision.This precision can meet the DEM cartographic chart of a scale 1:200 000, and in some certain key area, of a scale 1:50 000.

DISCUSSION&CONCLUSIONS
The errors in waterline heightening may be attributed to: errors in the delineation of the shorelines from various temporal data; errors in the georeferencing of the delineated shorelines; errors in the heights predicted by the hydrodynamic model; and further, the large error in this case might due to the delineation of the shorelines at extremely turbid mudflat, semi-range of tide.
However, there are still some issues needed further study, including using sensors with a large ground resolution.
Exploring algorithms tuned for high concentrations of various substances and the local specific optical properties of these substances.And a simultaneous detection of water color and land-water boundaries is needed.Also a short time lag between acquisition of remote sensing and in situ data used for validation is required, and the last but not the least, the sufficient geophysical and ecological knowledge of the area with the expert are needed.
The paper explores the potential for using the waterline method with multi-spectral images with middle spatial resolution to construct a Digital Elevation Model of a wide and grading mudflat.A new intelligent correlating model of waterline detection considering different tidal stages and local geographic conditions proved to detect waterlines more efficiently and exactly.The spatial structure of elevation superimposing on the points of waterlines was firstly constructed and a more accurate hydrodynamic ocean tide grid model was used.By the newly constructed abnormal hydrology evaluation model, a more reasonable and reliable set of waterline points was acquired to construct a smoother TIN and GRID DEM.This technique proved to be an effective tool to help coastal zone management and environmental protection.

S
is the scope of the tide gauge area t refers to time parameter In the equation (1), instant water level field 3 and Figure.4,the Co-pulse lag lines of M2 and S2 component tide (dashed) are not vertical to the coastline but have an oblique angle to the coastline.

Figure. 3
Figure.3 Co-amplitude line and Co-pulse lag line of M2 component tide

SS
=ordinal number ranked by the theoretical resolving model from high level to low level 2 =ordinal number ranked by the actual geographic and meteorological conditions from high level to low level N =the number of the waterlines n =the weight assignment value Results showed that the weight value had positive correlation with the in-situ surveying value (Figure.6).

Figure. 6
Figure.6 BJ-1 Waterline elevation abnormal scatter chart (a)Kriging interpolated DEM (b) TIN generated GRID DEM Figure.8Two methods of DEM generating Two DEMs of different intervals were built by this approach and were compared to evaluate the morph-dynamics of mudflats (Figure.8).TIN and Kriging interpolated methods have reasonable and expressive morph-dynamics distribution, but Kriging was better in details than TIN.From Kriging interpolated DEM, we continued to used ArcGIS Analysis Tool to generate the slope map and aspect map as follows Figure.9 and Figure.10.

Figure. 9
Figure.9 Slope of Mudflat (%)Figure.9showed the cartographic fruit of the slope of mudflat of Bohai Bay.The slope was expressed as percent or the ratio of elevation increment and horizontal distance.Results showed that the overall slope of the Bohai Bay tidal flat was gently sloping, between 0.016% and 0.08%, the average slope was 0.48‰.Analysis from the remote sensing images and the chart

Figure. 10
Figure.10 Aspect of T Mudflat (°) Figure.10 showed the cartographic fruit of the Aspect of mudflat of Bohai Bay.Aspect expresses the maximum variety direction of the elevation change.The North aspect is 0°and calculating in clockwise, the value scope is 0°~360°.In ArcGIS software, there are 9 aspects: tiny slope(-1), north slope(0°2 2.5°，337.5°~360°),north-east slope(22.5°~67.5°),east slope(67.5°~112.5 °), south-east slope (112.5°157.5 °), south slope(157.5°~202.5 °), south-west slope(202.5°~247.5 °), west slope(247.5°~292.5 °) and north-west slope(292.5°~337.5°).As for the results, lots of the aspects of the gradient were eastern and north-eastern.So the overall slope leant from land to sea and from south to north.The geographical explain was that the phenomena resulted in the tidal dynamical character and tidal flat aggradations character.Furthermore the aspects were caused by Qikou hydrology transect and Bohai Bay morphology sediment and erosion.The mechanism of the geography showed that the mudflat was mainly influenced by the sand transportation from the southern Huanghe River(Wang Y., 2012).So the aspect leant from south to north.This study not only gave a quantitative prove to the geographic discipline, but also filled up the mapping gap for this area.

Figure. 13
Figure.13 Histogram of error statistics of inversion and LiDAR DEM

Figure. 14
Figure.14 Inversion elevation transect compared with ground leveling with a horizontal interval of 60 m Compared the in-situ surveying data of the three transects and the inversion elevation data, the differences were partly high and partly low, the overall precision of the middle and high tide flat was better than the low tide flat.The variable trend was the descending from the land to the sea (Figure.14).

Table . 2
Statistic analysis of relative and absolute error between inversion and RTK surveying