Geophysical Mapping of Subsurface Stratigraphy Beneath A River Bed Using Ground Penetrating Radar: Lagos Nigeria Case Study

Ground penetrating radar (GPR) has been found to be most suitable for subsurface stratigraphic and lithologic mapping. In this study, we showed its application in bathymetric study in the region where a proposed pedestrian bridge across a river channels is planned. The application of a novel digital signal processing approach applied to the GPR acquired is shown. The methodology adopted involved the acquisition of GPR data along a 100 m transect, across the river using a 200 MHz GSSI SIR equipment. The GPR data was processed using the sequence viz: 1) Dewowing 2) Filter in frequency domain 3) Custom Gain 4) Automatic gain control 5) Delete Mean trace and 6) Deconvolution. In addition, complex processing methods (Instantaneous attribute) was applied to increase interpretability of the GPR data. The performance of each processing step is compared by examining its result. The study shows that filtering in frequency domain; gain application and delete mean trace are the basic processing steps used in low noisy data acquired with GSSI equipment. From the processing, it was evident that deconvolution is very effective in removing reverberation and/or multiples from the dataset. Delete mean trace produced the best result in removing high amplitude, continuous horizontal reflections (such as direct wave) that conceals true reflections. Even though these methods were effective, each method has its own advantages and limitations. Therefore, preservation of the geological features should supersede any other motive or decision made during GPR data processing. Among the complex processing methods, Instantaneous phase proved to be most resourceful for bathymetric study. The Interpreted section revealed the base of the river at a depth of 5 m. A good knowledge of the river geometry which will help civil engineers in proper planning of the bridge footings structure was also obtained.


Introduction
Rivers have a prominent role in many contexts as a natural environment, as a transfer medium, as a physical medium, as a natural resource. This list is not exhaustive. Therefore, in preserving rivers and ensuring quick transportation of man and goods, bridge construction becomes necessary. The process of designing a bridge however requires detailed knowledge of riverbed topography (geometry of the riverbed at various spatial and temporal scales).
Remote sensing technique such as Ground penetrating radar (GPR) and Bathymetric Lidar has been used to map river bathymetry [1][2][3]. Ground penetrating radar is a non invasive, non destructive electromagnetic (EM) investigation method. It is mostly used in reflection mode where a short pulse EM wave is emitted via an antenna into the subsurface. The arrival time and the amplitude of the reflected signals caused by changes in materials electrical property (Dielectric permittivity and Electrical conductivity) is recorded and analyzed. It provides information about both the spatial and vertical variation of the subsurface material at accuracy of few millimeters to several meters. Its high resolution power and ease in data acquisition has made it a preferred method. GPR can also be used in conjunction with other methods in bridge construction [4]. The use of high quality GPR data in mapping soil and rock stratigraphy are well documented in literature. For example, Travassos and Adepelumi [5] successfully carried out a combined marine and land GPR survey on the Jacui river in southern Brazil using a Pulse EKKO IV radar system with 50 MHz antennae and a 1000V pulser. Adepelumi and Olorunfemi [6] conducted an integrated engineering geological and geophysical investigation of the reclaimed Lekki Peninsula, Lagos, Southwest Nigeria where the subsurface stratigraphy in the eastern part of Lagos was delineated.
This study is aimed at helping civil engineers in the proper planning of bridge footing with understanding of the river geometry. The objective of the study includes; 1) 2D GPR data acquisition over the river; 2) to derive a processing step that will give a section that best represents the subsurface geometry; 3) mapping of the river geometry for bridge footing; 4) mapping of structural features present and 5) mapping of the water table.
Ground penetrating radar (GPR) is extensively used for a variety of applications. Among all the high-resolution geophysical methods, GPR has proven to be the most suitable for detection of karstic cavities, stratigraphic mapping [7] and bathymetric study, in a wide range of soil and rock conditions. However, one of the main limitations of GPR is the exact determination of the mean velocity of the electromagnetic waves, which is a key datum to estimate the depth of penetration into the ground.

Geology and Hydrogeology
The investigated area in located in the eastern part of Lagos State Nigeria. The area lies within the Dahomey basin; and it is a zone that is having coastal plateau geomorphic unit [8]. The area is within the rain forest zone of Nigeria with mean annual rainfall of above 1500mm. The climate of the area has two regimes of wet and dry season controlled by the dry north easterly air mass of Sahara origin and the humid maritime air mass blowing from over the Atlantic, both blowing from nearly opposite directions. The surface water network drains into the lagoon and Atlantic Ocean Southwards.
The Dahomey basin dates back to early Cretaceous after the separations of African -South America landmasses and subsequent opening of the Atlantic Ocean [9]. The oldest sedimentary unit on shore within the basin is the Cretaceous grits and sandstones of Abeokuta group. These strata are overlain on shore by the marine Paleocene limestone and shale of the Ewekoro formation. Overlying the Paleocene strata are the sand and clays of Ilaro formation. These lower Tertiary beds are overlain by the Coastal Plain sand composed of transitional to continental sands of Oligocene to Recent age.
The surface local geology of the study area is consistent with the regional geological setting of Lagos area as described by Longe et al. [10]. The surface geology of the study area is made up of the both the Ilaro Formation and the Recent littoral alluvial to Coastal Plain sand deposits (Benin Formation). The local subsurface geology reveals a lateritic cover, reddish brown in colour with sand and clay portions. The lateritic cover is of variable thickness. This lateritic cover overlies an alternating sequence of sand and clay deposit. Underlying this is a water-bearing zone consisting of loose, medium to coarse sand. The Ilaro Formation consists of rather massive sandstone with local clay intercalations. The Ilaro Formation is fine to medium grained, and is fairly well sorted. The The Ilaro Formation lies conformably on the Oshoshun Formation (Lower -Middle Eocene) and locally unconformably beneath the Benin Formation (Oligocene-Pleistocene). The Ilaro Formation is mostly likely to be Middle to Upper Eocene in age. The Benin Formation consists of continental sands with shale intercalations usually with good groundwater potential. The Ilaro Formation is estimated to be about 70 m thick and shows rapid lateral facies changes. This can affect the aquifer quality of the Ilaro Formation [11]. However, the underlying Ewekoro Formation is known to have good groundwater aquifer (limestone).
The water-bearing aquifers consist of sands, gravels or a mixture of the two. Textural variations from fine through medium to coarse sands and gravels occur and they are poorly to well sorted. When near the surface, the sand deposits are generally loose but become moderately dense with depth and occasionally with clay interbeds [10].

GPR Data Acquisition
For this study, GPR data was acquired across a river 62 m wide. The GPR data used was collected using a Geophysical Survey System Incorporated SIR system-3000 equipment. The survey was carried out using a 200MHz monostatic antenna with the antenna oriented parallel to the survey direction (Parallel-Broadside). The orientation of the traverses was approximately north-south with survey direction of 182ᵒ. Traces were recorded with a stack of 16 fold enough to improve the signal to noise ratio. 33 scans per meter were collected (0.03m station spacing) with a sampling window of 2000ns with an offset of +200ns. The GPR data positioning was calibrated using a survey wheel. Each radar traces contains 512 points per trace. The river line data was acquired with the antennae mounted on a wooden boat with the survey wheel manually rolled as the boat moves.

Results and Discussion
For proper interpretation of GPR data, the noise which masks the real reflections must be removed. The act of removing this noise is termed as GPR processing. Just as in seismic method, GPR processing is aimed at improving the interpretability of radar sections. This is done by enhancing the desired weak signals at the expense of the unwanted signals (noise). This increases the signal to noise ratio. This is however non unique and also depends on the objective of the survey. This section describes the processing flow applied to the acquired GPR data across a river.
Several authors have published different works on GPR processing using different steps and software in different environment. Fisher et al. [12] applied two separate but similar processing sequences for two different GPR data acquired using different antenna frequencies. Processing operations performed on the first data (using 100MHz frequency) include signal saturation, gain recovery, spiking deconvolution, bandpass filtering, and normal moveout corrections. A failure surface with no surface expression is 12 Geophysical Mapping of Subsurface Stratigraphy Beneath A River Bed Using Ground Penetrating Radar: Lagos Nigeria Case Study recognized in the final section. The second profile is a single-fold, 50 MHz, 0.5km line showing deltaic crossbed structures. These were processed using similar steps as above as well as f-k migration. Kim et al. [13] proposed novel techniques and sequence for removing ringing noise from GPR data. The proposed sequence by [13] includes background removal, f-k filtering, predictive deconvolution with filtering in wavenumber domain, and filtering by radon transform. Lee, in his quest to analyse the morphology of shoreface deposits in the Upper Cretaceous Ferron Sandstone Member, Utah, used an inhouse, well-established GPR processing system that includes 3-D depth migration [14]. Migration was performed with a typical dry sandstone velocity of 0.12 m/ns (0.39 ft/ns), producing a depth-migrated 3-D radar cube.
In this study, we have reviewed basic signal processing techniques and Instantaneous attributes. We showed the processing steps followed in removing the unwanted signals (noise) from the GPR data acquired over a river channel.

Dewowing
This is also known as signal-saturation correction. The short time interval between transmitter 'shots' during surveying and large energy input from the airwave, ground wave and near-surface reflections causes the receivers to become signal saturated [15]. Thereby causing a slowly decaying low-frequency 'wow' to be induced on the trace, that is, superimposed on higher frequency reflections. This can be corrected by increasing the antenna separation before acquisition or by introducing a high pass filter to cut the low frequency wow as it was done here. The effect of this correction when applied to the GPR section was not visible (Fig. 1)

Time Zero Shift Correction
In order to correct for misalignment of the first break that can occur in a radar reflection profile. GPR software allows the automatic realignment of GPR traces. They do this by physically moving individual or large sum of traces up or down by the required amount of two way time (TWT). As a threshold is required to identify the first break, problems can occur when data are particularly noisy, because it may be difficult to define consistently across all traces [16]. For this survey, a time shift (200 ns) applied during acquisition (to permit the visualization of the first arrival time) was the major target of our shift (Fig. 1). This process effectively removed all the unwanted section above the first arrival time and the airwave as shown in figure 2. This enables the interpreter calculate the true depth to subsurface reflector(s) as against whatever depth measured before this processing step.

Filter in Frequency Domain
Filters are designed to alter the shape of individual traces through mathematical manipulation, enhancing or eliminating certain features [17]. Filter is used when there is a need to enhance or reject a specified spectrum component of a signal. Therefore, this processing step is a basic requirement in every GPR processing sequence. Filtering in the frequency domain can be described by the equation, Y(ω ) = H(ω ) ⋅ X (ω ) (1) where, X (ω) is the spectrum of the raw signal, Y(ω ) is the spectrum of the processed signal and H(ω) is the spectrum of the filtering function.
This process is normally much faster than convolution, because we can use FFT in Fourier transformation. Band pass filter was applied to the section in figure 2. The frequency band chosen for filtering was 70/90 -150/220 MHz. Figure 3a-c shows the frequency spectrum obtained from unprocessed section, the obtained section after dewowing and the filtered section. Signal from the specified frequency spectrum was enhanced at the expense of others farther from the centre frequency, thereby giving a better image of the subsurface structure between the time windows 0 and 100 ns and 650 ns and 900 ns (Fig. 4). Observed on this section is the, air wave, river slope, few cultural features and noise. Although the wanted reflection has been enhanced by screening out unwanted frequency spectrum, certain processing steps are still needed to fully enhance the desired image.

Gain
Generally, radar signal strength reduces with time due to attenuation. Therefore, the possibility of having a strong reflection as the wave travels deeper into the subsurface is low. Amplitude of the weak signals in this part of the section can thus be increased by adding gain to the data. If the GPR recorded signal is f(t) and the amplifying function is defined as g(t) . The processed signal y(t) is given by: It should be noted that the multiplication is carried out between the time domain functions; therefore this is not a linear operation. The processing software in use, Radpro, has three different types of gains, which are automatic gain control (AGC), gain (also known as spherical exponential compensation, SEC) and custom gain. Where continuity of stratigraphic horizons is of central interest, it is often desirable to show all recorded information in traces, irrespective of amplitude. Consequently, in many sedimentological studies data are plotted and interpreted with automatic gain control [16].

Custom Gain
This is user specified. It is used to enhance an area of interest without affecting other areas with high reflections. This step was aimed at enhancing the time window between 30 ns and 80 ns on the radar section, making the weak subdued signals to be visible to the human eyes without over enhancing the other parts of the section. The enhanced section (Fig. 5) revealed the attenuated river bottom geometry between the time window of 20 and 64 ns.

AGC
When small reflection in the late time on a radar section cannot be clearly observed due to radar attenuation with time, AGC is applied. Such attenuation can be caused by factors including radiation pattern of antenna, reflection and refraction, and conductive loss. Therefore, AGC is often used to display GPR profile. Automatic gain control equalizes the energy spectrum along the whole section making it have uniform reflectivity power. It is mostly used when stratigraphic study is the main aim of the survey. However, we found it useful in revealing the upper part of the GPR section (Fig. 6). This confirms its resolution increasing power and its applicability in stratigraphic mapping.

Delete Mean Trace
One of the most common operations specifically applied to GPR is the use of average trace removal/ delete mean trace. This is very effective in allowing subtle weaker signals which are lost to become visible in a processed section. When the source antenna is placed on the surface, spherical waves are radiated both upward into the air and downward into the soil. Because of the continuity requirements for the electromagnetic field at the soil surface, the propagating spherical air wave gives rise to a lateral wave front in the soil. Similarly, the spherical wave propagating in the soil gives rise to the ground wave. The ground wave amplitude is known to decrease strongly with distance above the soil surface.
14 Geophysical Mapping of Subsurface Stratigraphy Beneath A River Bed Using Ground Penetrating Radar: Lagos Nigeria Case Study Raw GPR signal is normally composed of a strong directly coupled signal d(t, x) between the transmitter and the receiver and the continuing reflection signals r(t, x) . When the ground condition does not change a lot along a survey line, d(t, x) does not change, when the antenna separation is fixed. This processing step sums up all the traces in n turns, averages it and subtracts it from the raw GPR signal. Delete mean trace has the capacity of removing continuous reflection such as the air waves/ground wave (Fig. 7) thereby removing the low amplitude signals at the surface     It also works well in cases of multiple when the real signal is a horizontal continuous reflection as seen in the section below.Note: that this may create artifacts, particularly if there exist horizontally extended coherent regions.

Deconvolution
Deconvolution is an analytical process designed to remove the effect of a previous filtering operation [18]. In both seismic and radar, deconvolution attempts to remove filtering effects resulting from propagation of a source wavelet through a layered earth, and the recording system response [16].  Although Annan [19]believes that deconvolution is difficult to apply systematically to GPR profiles and often does not lead to a major improvement in resolution. However, it can be useful where reverberation is a significant problem. Predictive Deconvolution was applied to the section in figure 7 after specifying required values for parameters such as; Predictive lag (nsec) = 10.000, Operator length (nsec) = 400.000, and White noise content (%) = 5.000. Figures 8a-b shows the results obtained by using the correct and incorrect parameters. Figure 8b records the applicability of Deconvolution in radar processing especially in removing reverberation/multiples (40 ns, 15 m -20 m and 73 m -77 m), though this can be tricky. It can lead to various data artifacts as seen in figure 8b (77 m).

Advanced Processing
Just as in seismic attributes, GPR attributes have been utilized to emphasize geologic targets and to define critical rock and fluid properties. Geological interpretation of GPR data is commonly done by analyzing the radar amplitude, phase, and frequency in GPR map or cross section. Therefore, these three simple attributes -amplitude, phase and frequency remain the mainstay of geological interpretation in GPR attribute analysis [20, 21 / 19, 20 old].
Consider a complex trace z(t) (Fig. 9) whose wave patterns appear as a helix that spirals around the time axis. The projection of complex trace z(t) onto the real plane is the real trace x(t); the projection of z(t) onto the imaginary plane is the Hilbert transform trace y(t). Instantaneous attributes -amplitude a(t), phase Ф(t) and frequency ω(t)can be calculated from a complex trace using the following equations. Instantaneous amplitude a(t) = �x(t) 2 + y(t) 2 (2) Instantaneous phase Φ(t) = tan −1 ( Instantaneous amplitude attribute was generated using the Radpro software to evaluate the reflectivity and signal strength of the data. It helped in characterizing major reflectors which is indicative of the river base. This is observed in Figure 10 as region of high amplitude, with thick band of dark color. Instantaneous phase attribute (Fig. 11) enhanced the river base features more than the simple processed data and other complex processing steps. Phase changes dramatically to nearly 90º at the edges of the observed prominent anomaly. It thus provides further strong evidence for sharp transitions (edges). This shows the ability of Instantaneous Phase attribute to identify flat lying undulose reflectors, thin layers and curved reflectors.
Instantaneous frequency attribute (Fig. 12), just as Instantaneous phase has the ability to identify flat lying undulose reflectors and curved reflectors. It has the ability to separate reflections that arrive at similar time. When applied to the acquired GPR data, the continuity of the river base was enhanced and the reflectors were revealed irrespective of amplitude fidelity. The interpreted section can be seen in figure 13.

18
Geophysical Mapping of Subsurface Stratigraphy Beneath A River Bed Using Ground Penetrating Radar: Lagos Nigeria Case Study Figure 13 shows the interpreted GPR section. From this section, a distinct 2D image of the river geometry (in green colour) and the bathymetric water level (in orange colour) with depth ranging between 0 m and 0.5 m is revealed. Due to the limitation of monostatic antenna as regards CMP survey, velocity model used in the time -depth conversion was obtained from the velocity spectrum of the data. The river geometry observed from the interpreted section obtained from Instantaneous phase attribute can be described as follow: Width of river = 65 m Maximum depth = 5 m Maximum depth position = 36.5 m from the northern riverbank.
Other features observed are the sand ridge which lies between 29 m and 37 m at a depth of 3 m from the sea surface and the inferred old surface which extends northward from the ridge.

Conclusions
Ground penetrating radar as a 2D subsurface imaging tool is an efficient method in locating the best position for bridge footing in the river. Based on the analysis of the acquired data, it is concluded that GPR can be a useful, cost-effective tool for estimating water depths and identifying/mapping in-filled scour features in a clastic fluvial environment. However, the GPR tool will not provide quality images of sub-bottom strata if significant clay is present. The basic processing step applied to GPR depends on the quality of the data acquired, aim of the survey and the desired result of the interpreter. This study shows that filtering and gain (amplification) is essential in processing GPR data. Filtering in frequency domain requires the knowledge of the centre frequency of the survey which is usually between ± 50 MHz of the antenna frequency. Also, we have been able to establish the applicability of deconvolution in GPR data processing. These methods, however, cannot be routinely applied but may be a good approach when the GPR data acquisition is less erroneous. Complex processing methods (Instantaneous amplitude and phase) on the other hand is useful in enhancing the interpretability of GPR section. Considering the results obtained, we can conclude that the instantaneous phase attribute is most efficient in stratigraphic study, river base mapping and mapping of water table.