Contribution of Seismic and Acoustic Methods to Reservoir Model Building

The seismic reflection method has the advantage of providing a picture of the subsurface in three dimensions (3D) with a regular grid. In high resolution seismic surveys, the size of the grid cell is of the order of tens of meters for horizontal distances, and of several meters for vertical distances. The initial success of the seismic reflection method was solely linked to structural interpretation: it only recognized geometric shapes, regardless of the content. Stratigraphic interpretation permitted a deeper understanding of seismic data. To be validated, seismic data must be tied to wells. The tying and calibration of seismic data are carried out with the use of geophysical measurements in wells. Well measurements are seismic data such as vertical seismic profiles (VSP) with a resolution comparable to that of seismic reflection and log data with a vertical resolution ranging from tens of centimetres to several meters. Acoustic logs establish the link between seismic data (surface and wells) and other logs. The acoustic log has the vertical resolution of the other logs (electric, nuclear, etc) but obeys the same propagation laws as the seismic methods, although it operates in a different frequency bandwidth. The purpose of this chapter is to show the contribution of seismic and acoustic methods to reservoir model building. After a description of the geophysical methods (acoustic logging and reflection seismics), we show, with field examples, how the geophysical data have been recorded and processed to estimate porosity and hydraulic conductivity of the studied geological formations.

shear wave (only in fast formations), the fluid wave, two dispersive guided modes: the pseudo Rayleigh waves (only in fast formations) and the Stoneley waves.Acoustic logging allows the measurement of the propagation velocities and frequencies of the different waves which are recorded by an acoustic tool.Velocities and frequencies are computed from the picked arrival times of the waves (P-wave, S-wave and Stoneley wave).For a clean formation, if the matrix and fluid velocities are known, an acoustic porosity log can be computed from the acoustic P-wave velocity log.The analysis of the acoustic waves recorded simultaneously on both receivers of the acoustic tool is used to compute additional logs, defined as acoustic attributes, useful for the characterization of the formation, such as amplitude, shape index, wavelength and attenuation logs.Attempts to predict permeability from acoustic data have been made by several researchers.The historical focus has been on predicting permeability from P-wave velocity (Vp) and attenuation.The attenuation of a formation is dependent both on the clay content and the rock (clean part of the formation).Laboratory experiments (Morlier & Sarda, 1971) have shown that the attenuation of a clean formation can be expressed in terms of three structural parameters: porosity, permeability and specific surface.Both theoretical and experimental studies have identified the relation between acoustic attenuation and petrophysical parameters: with δ: attenuation (dB/cm) , f: frequency (Hz), ρ f : fluid density, φ: porosity , μ: fluid viscosity (centipoise), S: Specific surface (cm2/cm3), C: calibration coefficient, k: permeability (mD).
To compute permeability from eq 1, it is necessary to measure the P-wave attenuation and the porosity of the formation and calculate the effective specific surface of the formation.Ida L. Fabricius et al. ( 2007) have found that specific surface with respect to grain volume ( Sg) is apparently independent from porosity.In an attempt to remove the porosity effect on Vp/Vs and mimic a reflected φ vs log (Sg) trend, they propose to use the following relationship between porosity φ, Vp/Vs ( Vs: S-wave velocity) and Sg: (2) where it should be observed that Sg is multiplied by m to make Sg dimensionless.The a, b, c coefficients have to be computed to minimize the fluctuations of the specific surface Sg.In practice we try to predict permeability from P-wave attenuation, acoustic porosity, Vp to Vs ratio, and P-wave frequency.The log thus obtained is proportional to permeability and must be calibrated on core measurements or pumping tests.

Seismic reflection method
Today, most surveys are three-dimensional (3-D) designs where source and receiver points are distributed on a region of the earth's surface.This distribution typically forms a grid pattern, composed of source lines and receiver lines.When source and receiver lines are orthogonal, the 3-D survey design is called an orthogonal design.Design considerations for an orthogonal design include the spacing of the source and receiver lines, and the spacing of the source and receiver locations along the lines.The signal sent from any given source point is monitored by multiple receivers on different lines.The set of receivers that monitor a single source signal are often referred to as the template or active spread.
The classical approach to seismic processing can be summarized in two main steps.The first step includes pre-processing of the data and the application of static corrections.The purpose of pre-processing is to extract reflected waves from individual shots, by filtering out the parasitic events created by direct and refracted arrivals, surface waves, converted waves, multiples and noise.It is intended to improve resolution, compensate for amplitude losses related to propagation, and harmonize records by taking into account source efficiency variations and eventual disparities between receivers.Static corrections, that are specific to land seismic acquisition, are intended to compensate for the effects of the weathering zone (Wz) and topography.Static corrections can be computed by conventional methods such as the Plus-Minus method (Hagedoorn, 1959) or the refraction tomography method.The refraction tomography method, which is used to solve the travel time tomography, considers the velocity model provided by the Plus-Minus method to be a good solution for the global search and applies a simultaneous iterative reconstruction technique (SIRT) for the local search to produce an even more accurate result.The optimization scheme uses a misfit function defined as the least-square error between the computed and picked first break times.The velocity model is ray traced and updated until the misfit function value reaches a value previously defined by the modeler.The inversion technique (Mendes, 2009) gives the velocity distribution in the weathering zone (Wz) and is an effective method to compute the static corrections in 2D or 3D geometries.Records are then sorted in common mid-point gathers or constant offset sections.
If the data are sorted in common mid-point gathers, the second processing step is the conversion of the seismic data into a migrated seismic section after stack.This second step includes the determination of the velocity model, with the use of velocity analyses, the application of normal move-out (NMO) corrections, stacking and migration.If dip values are sufficiently large, velocity analyses provide time-depth relationships that are affected by dips.To overcome this inconvenience, a correction for dip effects or dip move-out (DMO) correction prior to determining velocity relationships must be applied to the data.If the data are sorted in constant offset sections, a pre-stack time (PSTM) or depth (PSDM) migration procedure which simultaneously performs dip correction, NMO correction, common midpoint stack and migration after stack, is applied.It is indispensable to have a good velocity model to carry out the migration process.The role of migration is to place events in their proper location and increase lateral resolution, in particular by collapsing diffraction hyperbolas at their apex.Proper migration requires the definition of a coherent velocity field, which must be a field of actual geological velocities in migrated positions.Several methods can be used to generate velocity models such as tomography or stereo-tomography methods.The migrated section can then be depth converted and transformed in pseudo velocity or acoustic impedance sections if well data (such as acoustic logs) are available.The procedures used to obtain acoustic impedance sections are often referred to as Model-based seismic inversions which require an a priori impedance model (obtained from well data) which is iteratively refined so as to give a synthetic seismic section to match the seismic section to be inverted.The final impedance model can be converted in porosity by using an empirical relationship between porosity and acoustic impedance established at well locations.

First field example: Acoustic logging and 3D seismic surveying in a carbonate formation
The 3D seismic data were recorded in France at the boundary of the Meuse and Haute-Marne departments in the vicinity of the Andra Center (National radioactive waste management Agency).The acoustic data were recorded in well "EST 431" located on the Ribeaucourt township, in the national forest of Montiers-sur-Saulx, 8 km North-North-West of the Andra Center.

Geological context
One of the drilling platforms, located in the center of the studied zone, was used to study formations ranging from the Oxfordian to the Trias.The analysis presented here concerns borehole "EST 431" and covers the Oxfordian formation.The objective of this borehole is to complement the geological and hydrogeological knowledge of this formation.This formation consists essentially of limestone deposited in a vast sedimentary platform (Ferry et al., 2007).The limestone facies, which vary from one borehole to another, are generally bio-detritic with reef constructions.In this formation, porosity ranges between 5 and 20% and "porous horizons" of a kilometric extension have been identified.As far as hydrogeology is concerned the observed water inflows are usually located in high porosity zones (Delay et al., 2007).The base of the Kimmeridgian shales was observed at -258.3 m (100 m ASL) and the base of the Oxfordian limestones at -544.3 m (-186 m ASL).During the drilling, water inflows were detected at -368m and -440 m.At the end of the drilling, the well was left in its natural water.

Acoustic logging
Figure 1 shows the 3 m constant offset section in the 333 -510 m depth interval, opposite the geological description.On the acoustic section, the refracted P-waves appear in the 0.6 -1.2 ms time interval, the converted refracted shear waves in the 1.2 -2 ms time interval, and the Stoneley wave in the 2 -2.4 ms time interval.On the acoustic section, we can differentiate: an event at 345 m showing a very strong attenuation of all the waves; an interval showing a very strong slowing down of the P and S waves (360 -375 m); a relatively homogeneous mid-level interval (375 -397 m); a level which stands out because of its strong variations in P, S and Stoneley velocities (397 -462 m); a very homogeneous zone below 495 m with easily identifiable P and S waves, and an image of alteration between 501 and 507 m. Figure 2 is a comparative display of acoustic logs (P-wave velocity, P-wave frequency, acoustic porosity), borehole diameter, gamma ray log and NMR porosity.A strong correlation between acoustic porosity and NMR porosity can be noted (correlation coefficient: 0.86).
The correlation coefficient between acoustic porosity and borehole diameter is high (0.84).The attenuation of the formation is computed from the amplitudes of the refracted Pwave acoustic signal in selected time windows.The energy ratio (between the two adjacent receivers of the acoustic tool) gives the attenuation (expressed in dB/m) of the formation.Two attenuation logs have been computed with two selected time windows.The first window is a short window centered on the first arch of the acoustic signal.The second is a large window centred on the first three arches of the acoustic signal.The correlation coefficient between the two attenuation logs is high (0.78).This allows the computation of an average attenuation log with its associated standard deviation which gives an indication of the uncertainty.The event at 345 m shows a strong attenuation, a high porosity value, a significant increase of the borehole diameter and a high gamma ray value.It has been interpreted as a shaly carbonate layer.P-wave attenuation and permeability are both strongly dependent on clay content.The gamma ray log was used to compute both shaliness and the shaliness corrected attenuation log in order to evaluate the formation permeability.
It was also used to compute an effective porosity log (Figure 2, bottom right).Figure 3 (left) shows the attenuation logs with their associated standard deviation before (a) and after (b) shaliness correction.Figure 3 (right) shows the Sg specific surface log and the acoustic permeability log calculated from eq. 1.The fluid viscosity μ and density ρ f have been assumed to be constant (μ = 1 centipoise, ρ f = 1 g/cm 3 ).The acoustic permeability log detects three permeable zones at 368 m, between 400 and 440 m, and 506 m.The permeable zone located at 506 m corresponds to a high value of conductivity and is characterized by a low porosity (6 %), a 10 dB/m attenuation, but a significant decrease of the P-wave frequency and of the specific surface.
During the drilling, water inflows were detected at -368 m and -440 m.At the end of the drilling, the well was left in its natural water.The hydraulic tests and conductivity measurements conducted later on did not confirm the inflow at 368 m seen during the drilling, but they validated the 400 -440 m and 506 m permeable zones detected by the acoustic logging.The absence of active hydraulic fracture in this borehole was confirmed by hydrogeological tests and acoustic logs.
The 3D seismic survey was performed to predict the porosity of a carbonate formation and to evaluate the hydraulic conductivity between the porous bodies of the aquifer.

Seismic acquisition and processing
The 3-D survey design is an orthogonal design.The active spread is composed of 8 receiver lines with 74 receivers per line.The distance between 2 adjacent receiver lines is 100 m.The distance between 2 receivers is 25 m.A single source line is fired per active spread.The distance between two source points is 25 m.The source is a vibroseis source generating a signal in the 14-140 Hz frequency bandwidth.The distance between 2 adjacent source lines is 100 m.The cell or bin size is 12.5 m x 12.5 m.The nominal fold is 37.The size of the area covered by the 3D is 4 km 2 .A conventional seismic sequence was applied to the data set.It includes amplitude recovery, deconvolution and wave separation, static corrections, velocity analysis, CMP-stacking and time migration.After migration, a model-based stratigraphic inversion (a priori impedance model obtained from well data) provides a 3D impedance model cube.Figure 4 is an example of an impedance section extracted from the 3D block.It also displays the porous layers (Hp1 to Hp4) associated with the Oxfordian carbonate.

From 3D impedance to 3D porosity
At well locations, the acoustic impedance log was calculated from the density and velocity logs.The porosity vs. impedance cross plot, displayed in figure 5 (top), was used to define a linear law between the two.The porosity is expressed in % and the acoustic impedance in (m/s).(g/cm 3 ).The cross plot was obtained by using density, acoustic velocity, and porosity (NMR) logs recorded in 6 wells.This empirical relationship between porosity and acoustic impedance was used to convert the 3D impedance into porosity within the Oxfordian layers (Figure 5, bottom).The porous layers are clearly visible.The 3D cube makes it possible to provide 3D imaging of the connectivity of the porous bodies.Core analysis was carried out to establish porosity vs. permeability laws and it demonstrated that the permeable bodies have porosity larger than 18 %.Consequently, a porosity cut-off of 18% was used to show the connectivity of the porous bodies (figure 6).A porosity cut-off of 21% was also applied to extract the porous bodies having the best hydraulic conductivity within the porous layers.

Second field example: Case study of a near surface heterogeneous aquifer
Many underground aquifers were developed as experimental sites during the past decade.These sites are designed for in-situ measurements and calibration of flow, transport and/or reactions in underground reservoirs that are heterogeneous by nature.The University of Poitiers (France) had a Hydrogeological Experimental Site (HES) built near the Campus for the sole purpose of providing facilities to develop long-term monitoring and experiments for a better understanding of flow and transfers in fractured rocks (Bernard et al., 2006;Kaczmaryk & Delay, 2007a,b ;Bourbiaux et al., 2007).

Geological context
The aquifer studied, 20 to 130 m in depth, consists of tight karstic carbonates of Middle Jurassic age.It lies on the borderline, named the "Poitou threshold", between the Paris and the Aquitaine sedimentary basins (Figure 7).The Hydrogeological Experimental Site (HES) covers an area of 12 hectares over which 35 wells were drilled to a depth of 120 m (Figure 7).Hydrogeological investigations show that maximum pumping rates vary from well to well and range from 5 to 150 m 3 /h.The top of the reservoir was initially flat and horizontal, 150 millions years ago, but has been eroded and weathered since, during Cretaceous and Tertiary ages.It is shaped today as hollows and bumps with a magnitude reaching up to 20 m.The building phase started in 2002 and up to now, 35 wells have been bored over the whole thickness of the reservoir.Most wells dispose of documented drilling records and logs of various nature among which are gamma-ray, temperature, and acoustic logs.In addition, two wells were entirely cored.To sum up, the aquifer responds fairly evenly to the hydraulic stress of a pumped well, with pressure depletion curves merged together in time and amplitude whatever the distance from the pumped well.This merging is assumed to be the consequence of a local karstic flow in open conduits that have been unclogged by the drilling and pumping works.The presence of pervasive karstic drains is supported by recent logs in the wells using optic imaging.Almost all the wells have shown caves and conduits cut by the walls of the boreholes with sometimes mean apertures of 0.2 -0.5 m.These conduits are mostly enclosed in three thin horizontal layers at a depth of 35, 88 and 110 m.Of course, these layers are intercepted by vertical wells and this potentially results in a good connection between wells and karstic drains.This connection is mainly controlled by the degree to which drains are re-opened in the vicinity of the well.In the end, it was considered crucial to better image the geometry of the reservoir with a resolution compatible with, on the one hand, the scale of a well and, on the other hand, the scale of the entire experimental disposal.High resolution geophysical tools seem well designed to undertake that kind of investigation.

3-D seismic acquisition and processing
Due to the limitations of the area, the length of the seismic line could not exceed 250 m in the in-line direction.In the cross-line direction, the extension of the area does not exceed 300 m.As a result, 20 receiver lines were implemented, with a 15 m distance between adjacent lines.
Figure 8 shows the map locating the seismic lines and the wells.In the acquisition of data a 48 channel recorder was used.An explosive source (25 g) was detonated and a single geophone (10 Hz) per trace was deployed.Such a source makes it easy to identify and pick first arrivals.A 5 m distance between two adjacent geophones was selected to avoid any spatial aliasing.A direct shot and a reverse shot were recorded per receiver line.Figure 8 shows an example of an in-line direct and reverse shot gather.Three shot points in the crossline direction were fired at distances of 40, 50 and 60 m from the receiver line under consideration.Figure 8 shows an example of a cross-line shot gather.The range of offsets was selected to optimize the quality of the seismic image over the reservoir depth interval, between 40 and 130 m.The minimum offset distance was chosen equal to 40 m to reduce the influence of the surface waves.The time sampling interval is 0.25 ms and the recording length is 0.5 s.

Seismic refraction survey
To obtain the velocity of the refractor (top of the reservoir) and its depth, the Plus-Minus method was used.This method requires recordings where geophones are aligned with shot points.The arrival times of the direct and refracted waves were picked on all the in-line shots.The procedure was applied on each line independently.In order to obtain a map with a homogeneous sampling interval both in cross-line and in-line directions, the delay time curves were interpolated by kriging (Chilès & Delfiner, 1999) on a grid 2.5 m x 5m.
To perform the depth conversion, the velocity of the medium (Wz) situated above the refractor must be known.Here, it is given by the slope of the direct wave.In the area, the velocity V2 of the refractor was found to be 3350 m/s and the velocity of the Wz to be 850 m/s.On the Wz depth map (Figure 8), the arrow indicates the direction N 90° which corresponds to the main orientation of fracture corridors.
In order to add information in the inversion procedure, we used simultaneously in-line and cross-line cross shots with an offset of 60 m.The shots were selected to be sure to have the refracted wave as first arrival wave whatever the source receiver distance.The picked times of the first seismic arrivals on all the shots (in-line and cross-lines shots), the Wz depth map and the velocity model obtained by the Plus-Minus method are input data for the inversion procedure.The 3D a priori model is a homogeneous two layer model of V1 = 850 m/s and V2 = 3 350 m/s with the interface given by the Wz depth map.
For the inversion, we worked with 56 shots simultaneously: 40 shots acquired along the 2Dlines plus 16 shots acquired with a cross-offset equal to 60 m.A 400 Hz frequency was considered to compute the Fresnel volume.This value was estimated after tests with several frequencies to study the sensitivity of model velocity to the smoothing nature of the Fresnel volume method.
As we consider the 3D a priori model a valuable approach for the background velocity, we used it for the inversion SIRT procedure.Figure 9 shows the result of this process which was stopped when the rms error reached a value around σ = 0.87 ms, with an improvement of about 80% during the optimization procedure.Figure 9 shows the velocity distribution at 20 m in depth and the 2500 m/s iso-velocity depth map.We can notice a strong correlation between the 2500 m/s iso-velocity depth map (Figure 9, left) and the Wz depth map (Figure 8, bottom left).The correlation coefficient reaches 0.96. Figure 10 shows a 3D block with velocity sections located at a 0 m, 60 m and 180 m distance in the cross-line direction and velocity map located at 20 m in depth.
The velocity distribution versus depth can be represented by a two-layer model.The interface corresponds to the top of the karstic reservoir.Above the interface, the velocities are low (Vp ranging from about 850 m/s to about 2500 m/s) with no significant lateral variation but with a strong increase in depth associated with the 2500 m/s iso-velocity in the vicinity of the interface.
The inversion results obtained with 3D data emphasize the geological structures mentioned previously in Mari et al. (2009) allowing a better recognition of their alignments and shape (corridor of fractures).Furthermore, no cavities were detected near the surface.The velocity model obtained by inversion of the first arrival picked times was used to extend the 3D velocity block obtained by the reflection survey (Mari et al. 2009) to the 0-35 m depth interval, emphasizing the necessity of employing 3D acquisitions for near surface studies.

Seismic reflection surveying
The processing sequence has been described in detail by Mari & Porel (2007).Each shot point (both in the cross-line direction and in the in-line direction) was processed independently to obtain a single-fold section with a sampling interval of 2.5 m (half the distance between 2 adjacent geophones) in the in-line direction and with a sampling interval of 5 m in the cross-line direction.For seismic imaging based on reflected waves, it is necessary to be able to separate weak reflected events from high-energy surface waves such as pseudo -Rayleigh waves.Wave separation is a crucial step in the processing sequence.This specific field case illustrates the benefit of combining two different wave-separation methods in order to remove all the energetic wave-field.The conventional F-K method was used to filter surface waves and converted refracted waves.The SVD method (Singular Value Decomposition) was then used to extract refracted waves.After amplitude recovery, deconvolution and wave separation, normal move-out corrections were applied to the residual sections in order to obtain single-fold zero offset sections at normal incidence.In well C1 situated in the central part of the site, a Vertical Seismic Profile (VSP) was recorded.VSP data were processed to obtain a time versus depth relationship and a velocity model.The velocity model was used to apply the normal moveout corrections.The VSP time versus depth law was also used to convert the time sections into depth sections with a 0.5 m depth sampling interval.The single-fold depth sections were merged to create the 3D block.The width of the block in the in-line direction equals 240 m and 300 m in the cross-line direction.In the in-line direction, the abscissa zero indicates the location of the source line.The abscissa of the reflecting points varies between -120 m and 120 m in the in-line direction, the distance between two reflecting points equals 2.5 m.In the cross-line direction, the distance between two reflecting points equals 5 m.The depth sections were deconvolved in order to increase the vertical resolution.They were then integrated with respect to depth to transform an amplitude block into a 3D pseudo velocity block in depth, using velocity functions (acoustic logs recorded at wells C1, MP5, MP6, M08, M09) as constraints.The pseudo velocity sections of the 3D block thus obtained were merged with those obtained by refraction tomography (Figure 10) to create a 3D extended velocity model from the surface.Figure 11 shows the in-line pseudo velocity sections extracted from the 3D extended velocity model and situated at cross-line distances of 0, 60 and 180 m.It also shows the velocity map at 87 m in depth.The 3D velocity model shows the large heterogeneity of the aquifer reservoir in the horizontal and vertical planes.

3D porous model building
In the area covered by the 3D seismic surveying, 23 wells (location shown in Figure 7) have been drilled.The wells are regularly spaced (∼ 50 m) and are used to perform many hydraulic tests (interference pumping and slugs).In the wells, several logs have been recorded (electrical and gamma ray logs).Due to the homogeneous spatial distribution of wells, in which resistivity logs have been recorded, we chose to use a method based on electrical measurements to quantify the 3D porosity distribution within that aquifer.The seismic interval velocity-to-porosity conversion was performed in two steps (Mari et al., 2009): • first step: from 3D interval seismic velocity to 3D resistivity • second step: from 3D resistivity to 3D porosity.Faust (1953) established an empirical relationship between seismic velocity V, depth Z, and electrical resistivity measurements Rt.For a formation of a given lithology, the velocity V can be written as a function of the depth Z and resistivity Rt as follows:  exponent), of that empirical law were determined through a least-square minimization of the difference between the 3D-block-extracted seismic velocities and the velocities predicted from Faust's law using the long normal resistivity log data as input.2D distribution maps of the C and b values over the site could then be built from the calibrated values in each of the 22 wells.These maps were used for the velocity -resistivity conversion of the 3D seismic block.Figure 12 shows the results obtained at well MP6 and the resistivity map at 87 m in depth.Archie (1942) showed empirically that for water-saturated permeable formations, the relation between the true formation resistivity, Rt, and the resistivity, Rw, of the water impregnating the formation was given by :

From 3D resistivity to 3D porosity
where F is the " resistivity formation factor ". φ is proportional to the formation porosity and m is a " cementation factor", that is a formation characteristic.The F value derived from the resistivity measurement, Rt, is unaffected by the mineralogical constituents of the formation matrix.Although the "cementation factor" value may vary between 1.3 and 3 according to the formation lithology, an approximate value of 2 is generally adopted for a well-cemented sedimentary log.Although applicability of Archie's law may be argued a priori for a karstic reservoir, two reasons motivated its choice.Firstly, the reservoir remains essentially a sedimentary carbonate formation at the seismic resolution scale.Actually, the size of the seismic bin (2.5 m in the in-line direction, 5 m in the cross-line direction ), and the seismic vertical resolution ranging between 1 and 2 m, lead to an elementary seismic cell volume of 12 m 3 at least.Secondly, the volume of the karstic bodies represents only a small percentage (2 to 3 %) of the reservoir volume.This volume was estimated by analyzing borehole images (Audouin & Bodin, 2007).For the above two reasons, the previous seismic-derived 3D resistivity block (Rt-seis) was converted into a 3D pseudo-porosity block, by using the Archie-law given by eq. 5 with m= 2 and with the resistivity of the formation water, Rw, estimated at 20 ohm.m. Figure 12 shows the porosity log computed from the seismic-derived resistivity log, Rt-seis, at well MP6.The production profiles measured in many wells revealed that these layers actually correspond to the major water feeding levels of wellbores.In order to further analyze the spatial distribution of porous bodies and of presumably-conductive flow paths, different cut-off values were applied to the 3D seismic porosity block.Figure 13 shows several 3D seismic pseudo-porosity blocks, associated with porosity cut-off values of 10 and 30%.The blocks were extracted from the central part of the area.They have a length of 120 m in the in-line direction and a length of 300 m in the cross-line direction.The extracted 3D reservoir volume with a porosity smaller than 10% (Figure 13, top right) actually represents the largest fraction of that aquifer, i.e. tight carbonates with a low permeability (less than a millidarcy).This is consistent with the observation of very sparse and channeled flow paths within that aquifer.We may indeed assume that, within that extracted low-porosity volume, the density of conductive (karstic) bodies is too low to ensure a hydraulic communication between wells, because the velocity-to-porosity converted block is derived from a high-resolution 3D seismic block, defined at a metric to pluri-metric scale that is significantly less than the well spacing.
The bodies with porosity higher than 10 % (Figure 13, bottom left) are mostly distributed within 3 layers, located at 30-35 m, 85-87 m and 110-115 m depth.Finally, a cut-off value of 30% was applied to the seismic porosity block to show up the most porous bodies of that aquifer.These highly-porous bodies (Figure 13, bottom right), mainly located in the intermediate porous layer situated in the 85 -87 m depth interval, represent only 2 % of the total volume of aquifer investigated in this study.

3D porosity and hydraulic conductivity
The hydraulic slug tests, performed at the HES, show a very rapid propagation of the pressure wave over large distances, say 100 m on average.These observations allowed us to map a diffusivity distribution and the importance of connections between wells (Figure 14).Preferential connections are visible along the N90 direction (wells M13-M21-M22-M19 and wells M04-M06-M11).Incidentally, logging data show three water producing layers (35, 87, 115 m depth, respectively).Their presence is not systematic however at each well of the site.
The 115 m-depth layer is located on the SW and NW borders of the site whereas the 87 mlayer is present everywhere.The variations of the seismic pseudo-porosity at different depths can be compared directly with hydrogeological data.The porosity map (Figure 14) drawn for the 87 m horizon confirms both hydrogeological data (flowmeters) and optical imaging obtained from the majority of wells.Actually, there is a significant connection between M13 and M21 which is also visible in the form of a high seismic porosity zone.It is therefore very likely that high porosity zones of geophysical maps do correspond to water productive areas.Interference testing gives information about hydrodynamic behaviour at the site scale but does not image the local hydraulic connection between wells, even though pressure transients may differ from one observed well to the other (Bernard et al., 2006;Kaczmaryk & Delay, 2007a).
It can be now wondered on the interest of a 3D porosity block for modelling flow in the aquifer.As stated above, the insurgence of karstic features is visible both onto the block and on responses to hydraulic stresses showing very rapid propagation of pressure wave depletions between distant wells.The question raised by these observations is how to model the static geometry of the porosity distribution and the dynamic hydraulic responses within a flow simulator that would keep inversion feasible.Stated otherwise, the point is to assign hydraulic parameters to the block by matching up numerical simulations onto flow data generally taken from hydraulic tests.In view of the spatial structure depicted by the porosity block, a continuous approach to flow seems well suited to fit a block of hydraulic conductivity (and specific storage capacity) that would satisfy the continuity of spatial porosity distribution drawn from seismic data.This makes us overlook discrete fracture networks and other object-based representations of fractured karstic aquifers, even though these interpretative depictions could probably be applied to the present case study.Most of the hydraulic tests performed at the HES were based on the notion of interference between wells.It is recalled that interference testing monitors and interprets the system's response in terms of hydraulic head (or drawdowns) variations in an observation well located at some distance from another well where a pumping (injection) stress is applied.
Observations are usually representing a state of the system integrated over the whole screened intervals of the wells.Several variations of this basic configuration are then implemented.www.intechopen.com These include, for instance, dipole-tests, where groundwater is circulated between a pumping and an injection point, or cross-borehole flow logging, where the observation well is monitored at several depths.In general, the collected information is therefore partly integrated along the vertical direction but relevant to continuous approaches to flow.It can be sometimes dwelled on the capability of interference data to really conceal all the elements necessary for the assessment of 3D flow features (Delay et al., 2011), but this technical discussion is beyond the scope of the present contribution.
A classical dual-continuum approach was recently modified by Kaczmaryk & Delay (2007b) to account for the rapid propagation of pressure waves in a fractured karstic aquifer.This approach separates two overlapping continua, the fracture medium and the matrix medium, respectively, and adds hyperbolic wave propagation to the fracture medium.The local flow equations are written as is the vector of spatial coordinates, h f (x, t) and h m (x, t) are hydraulic head in the fracture and in the matrix continua, respectively; K f (x) [LT -1 ] represents the hydraulic conductivity of the fracture medium, Ss f (x) and Ss m (x) [L -1 ] respectively represent specific storage capacity in the fractures and in the matrix, α(x) [L -1 T -1 ] is the exchange rate between fractures and matrix.From eq. 5, three basic types of models can be analyzed: (1) a single medium approach, i.e., u = 0, Ss m = 0 and α = 0; (2) a classical dual porosity medium, i.e., u = 0 and Ss m ≠ 0, α ≠ 0; and (3) a modified dual porosity medium, i.e., u ≠ 0 and Ss m ≠ 0, α≠ 0. These three types of models can be either taken homogeneous or heterogeneous in space.The first assumption is well suited to interpret drawdown curves one at a time (i.e., the single relationship between the stress at the pumped well and the response at the observed well).By repeating inversion runs of a homogeneous model for all pairs pumped -observed wells available at the HES, one obtains a set of hydraulic parameters that bound the reasonable (or acceptable) values for the bulk hydraulic behavior of the system.The geometric structure of the aquifer revealed by the porosity block is obviously overlooked but this homogeneous approach may fix the ideas for further modeling tasks in the heterogeneous context.For instance, it could be considered that i-high porosity bodies are highly permeable with hydraulic conductivity values on the order of the upper bound of the homogeneous approach, ii-low porosity areas are affected with the small values of hydraulic conductivity from the lower bound of the homogeneous approach.The exercise of inverting all drawdowns curves one at the time was carried out by Delay et al. (2007), Kaczmaryk & Delay (2007a, b) for the three main homogeneous configurations of eq. 5, i.e., a single continuum, a dual fracture-matrix continuum, a dual continuum with wave propagation.In the end, it was shown by the authors that the relationships between local stresses and local responses were sensitive to the three assumptions prevailing to flow but always rendered quite homogeneous hydraulic conductivity values.Actually, the three homogeneous configurations differ by their capacity parameters (specific storage and, to some extent, the exchange rate between fractures and matrix) that may vary in any case over wide ranges.Another attempt of drawdown inversion was carried out by Ackerer & Delay (2010).The authors considered flow in a single continuum model averaged over the vertical direction, thus yielding a diffusive 2D model widely used in groundwater engineering.They also inverted fictitious scenarios adding sequentially data obtained during various interference tests with the aim of providing to the calibration procedure several stress-response pairs evenly spread over the system.A result in terms of drawdown fitting and sought map of hydraulic conductivity is reported in Fig. 15.Again, there is no prior conditioning onto the structure given by the porosity block but the inverted hydraulic conductivity field is structured in space, obeying an exponential covariance not far from that used to filter seismic data.As shown in Fig. 15, inverting several drawdown curves results in variations in space of the hydraulic conductivity over four orders of magnitude, 10 -8 -10 -4 ms -1 .This range could be considered as a good forecast of values that could be assigned to a 3D conductivity block inheriting from the structure provided by seismic data.
The exercise of inverting a 3D block of hydraulic conductivity conditioned onto the geometry of the porosity block is under investigations.The targeted objective is giving to the inverted fields of hydraulic parameters some consistency with the 3D spatial structure of the main flowing bodies.Several options can be proposed to deal with a tractable problem.The first option, which is quite inescapable, is that of thresholding the sought values of hydraulic conductivity.To make it short, it can be proposed to separate the system in e.g., 3 types of conductivity values: high, mid and low values corresponding to similar cut-offs in the porosity values of the block.Without this simplification, it is likely that inverting a continuous distribution of hydraulic conductivity would become cumbersome or without tractable solution because of the too large number of freedom degrees given to the calibration procedure compared with available flow data.In the same vein, it can be discussed on the interest of proposing models of the dual-continuum type, eventually amended by adding wave propagations, in the case we dispose of a good image of the structural heterogeneity of the medium.In other words, is the picture given by the porosity block separating main flowing bodies from less permeable areas enough to represent the complexity of flow in a karstified aquifer, or is it needed to include locally more physics in the flow equations?Tackling the problem needs simply increasing the complexity of the local flow equations solved in a 3D block respecting the structure and then proceed with a sensitivity analysis of the models to parameters.The enhancement : single continuum → dual continuum → dual-continuum + hyperbolic wave propagation (see above) would inform on the relevance of adding or leaving out some mechanisms in the flow equations for a better understanding of measured data.Conversely, the sensitivity analysis while degrading the model would report on the capability of data to conceal information on the different mechanisms ruling flow in the fractured karstic aquifer.Notably, this topic remains a key question, unresolved for the moment, because very often the same set of data can be inverted with the same accuracy for different conceptualization of flow in the system.
Incidentally another problem hampers for the moment the passage from a 3D block of seismic velocities to its equivalent in terms of hydraulic conductivity distributions.As stated above it is quite straightforward to build a prior estimate of the conductivity block by copying it on the distribution of seismic velocities.But it is not sure at all that this block would result in a reasonable fitting of interference data.A close look at interference tests performed over the HES shows that they mix reciprocal and non reciprocal responses.It is reminded that reciprocity, as defined by Lorenz with reference to electromagnetic fields, states that for a given process in a heterogeneous medium (e.g., Darcian flow in a porous medium) a stress at a location A yields a response at location B which coincides with the response A due to an equivalent stress in B. In the present case of the HES, pumping a constant flow rate at well A may generate transient pressure drawdowns in well B that differ from the pressure response in A due to pumping in B. In these conditions, a 3D flow simulator for a single continuum based on the mass balance principle and Darcy's law should be discarded since it can be proved that this physical configuration always yields reciprocity irrespective the degree of heterogeneity in the medium and the boundary conditions (Delay et al., 2011).For 2D applications as discussed above, reciprocity gaps are damped by the integration of flow along the vertical direction making that hydraulic conductivity is integrated over the wetted thickness of the aquifer.For 3D flow with the aim of seeking accurately whether or not high-resolution seismic data coincide with hydraulic conductivity distribution, using non-reciprocal data for inverting a model which renders reciprocity would probably yield artifacts hardly graspable at this stage of our knowledge on the problem.Current studies (Delay et al., 2011) are analyzing the possible source of non reciprocity in the existence of inertial effects, non linear behavior as regard flow of the karstic features riddling the aquifer, complex flow as that in dual media, etc.These peculiar behaviors of the HES aquifer are incentive to first revisit the physics of a 3D flow model before proposing a block of hydraulic conductivity conditioned on both high-resolution seismic data and hydraulic tests.

Conclusion
We have shown with field examples how seismic and acoustic methods can be used to estimate porosity and condition inversion of hydraulic conductivity fields of geological formations for 3D models of reservoir.
The first field example concerns the study of t h e p o r o u s l a y e r s a s s o c i a t e d w i t h t h e Oxfordian formation.A high resolution 3D seismic survey covering an area of about 4 km² was recorded.The processing gave rise to a time migrated 3D block.After migration, a model-based stratigraphic inversion provided a 3D cube of impedance after calibration at several wells.A porosity vs. impedance relationship, obtained by using density, acoustic velocity, and porosity ( NMR) logs recorded in the wells, was used to convert the 3D impedance into porosity.The 3D cube makes it possible to provide 3D imaging of the connectivity of the porous bodies.A porosity cut-off of 21% was applied to extract the porous bodies with the best hydraulic conductivity within the porous layers.In order to complement the hydrogeological knowledge of the Oxfordian formation and to evaluate the power of the acoustic method to predict permeable zones, a full waveform acoustic log was recorded in an experimental well.The permeable zones detected by the acoustic logging were validated by hydraulic tests and conductivity measurements conducted later on .The absence of active hydraulic fracture in this borehole was confirmed by hydrogeological tests and acoustic logs.
The second field example concerns a near surface karstic reservoir.Both surface seismic refraction and seismic reflection data were recorded on the experimental hydrogeological site that has been developed for several years near Poitiers.Seismic reflection data enabled us to derive a 3D seismic pseudo-velocity block.Refraction tomography applied to the field data is conducive: (i) to obtain a complete velocity model in the first 35m; (ii) to map the top of the karstic reservoir; (iii) to detect the main corridor of fractures; (iv) to be sure that no cavities are present in the first 35 m.The velocity model obtained by refraction tomographic inversion was used to extend the 3D velocity block obtained by reflection survey to the 0-35 m depth interval.The 3D seismic pseudo-velocity block reveals a large heterogeneity of the aquifer reservoir in the horizontal and vertical planes.The low-velocity areas were found to correspond to the conductive levels and regions, as identified from well logging and flow interference tests.In order to quantify the porosity variations within that aquifer, the seismic-interval velocities were first converted into resistivity values.For that purpose, the empirical relationship between seismic velocity and true formation resistivity proposed by Faust (1953) was used.Resistivity values were then converted into porosity values, by using Archie's law (1942).The resulting 3D seismic pseudo-porosity block revealed three highporosity, presumably-water-productive, layers, at depths of 35-40, 85-87 and 110-115 m.The 85-87 m-depth layer is the most porous one, with bodies of a porosity higher than 30%, that represent the karstic part of the reservoir.That seismic pseudo-porosity distribution appears to be consistent with the available hydrogeological data recorded at the site.It is therefore very likely that the high seismic-porosity zones (higher than 30 %) of the geophysical 3D block correspond to water producing layers.However, the conversion of that pseudoporosity block into a dynamic-flow-property block is under investigations.In the case of the HES showing evidences of karstic features yielding complex flow and sometimes nonreciprocal behavior, the passage from seismic velocities to hydraulic conductivity is less straightforward as simply copying the block of conductivity on the block of velocities.However, we point out that high resolution seismic profiles inform with precision on the location of, e.g., conduit flow compared with matrix flow.To conclude, the very high resolution seismic survey of that near-surface aquifer makes it possible the construction of a 3D seismic-porosity block providing a new insight to a deterministic high-resolution model of reservoir.

Acknowledgment
We thank Andra for permission to use the data presented in the first field example.We thank Daniel Guillemot and Beatrice Yven for their valuable help and advice.

Fig. 1 .
Fig. 1.Acoustic logging and lithological description (Courtesy of Andra)The a, b, c coefficients of equation 2, computed in the 370 -500 m depth interval (a = 0.02, b = 0.012 and c = 6.25), were used to calculate the specific surface Sg with respect to grain volume.Figure3(right) shows the Sg specific surface log and the acoustic permeability log calculated from eq. 1.The fluid viscosity μ and density ρ f have been assumed to be constant (μ = 1 centipoise, ρ f = 1 g/cm 3 ).The acoustic permeability log detects three permeable zones at 368 m, between 400 and 440 m, and 506 m.The permeable zone located at 506 m corresponds to a high value of conductivity and is characterized by a low porosity (6 %), a 10 dB/m attenuation, but a significant decrease of the P-wave frequency and of the specific surface.During the drilling, water inflows were detected at -368 m and -440 m.At the end of the drilling, the well was left in its natural water.The hydraulic tests and conductivity measurements conducted later on did not confirm the inflow at 368 m seen during the drilling, but they validated the 400 -440 m and 506 m permeable zones detected by the acoustic logging.The absence of active hydraulic fracture in this borehole was confirmed by hydrogeological tests and acoustic logs.

Fig. 3 .
Fig. 3. Permeability estimation from acoustic logs (Specific surface, porosity, P-wave frequency).Left: Attenuation logs and their associated Std before (a) and after (b) shaliness correction, Right: Specific surface log and Predicted Permeability from acoustic logs.

Fig. 7 .
Fig. 7. Hydrogeological Experimental Site in Poitiers: site map and well locations.
top left: Seismic line implementation and well location (red points), top right: Example of direct and reverse in-line shot points, bottom left: Wz depth map, bottom right: example of cross-line shot point.

Fig. 8 .
Fig. 8. Seismic acquisition: examples of shot points and Wz-depth map obtained by seismic refraction surveying.
3)with V the P-wave velocity of the formation in m/s, Z the depth in m, Rt the electrical resistivity in ohm.m,C and b the coefficients associated with the Faust's equation.At each well where a long normal log was recorded, an interval velocity log was extracted from the 3D seismic interval velocity block.The two sets of data (resistivity and seismic velocity) were combined to calibrate an empirical Faust's law, which was then used as a local constraining function to transform the 3D pseudo-velocity block into a 3D pseudoresistivity.For each well, the two coefficients, C (constant coefficient) and b velocity sections located at a 0 m (a), 60 m (b) and 180 m (c) distance in the cross-line direction and velocity map located at 87 m in depth (d).
maps located at 87 m in depth (a, d).Velocity, resistivity and porosity logs at well MP6 (b).Velocity and porosity sections located at a 180 m distance in the cross-line direction (c).

Figure 12
Figure12also shows the in-line pseudo velocity and pseudo porosity sections (extracted from the 3D velocity model and from the 3D porosity block) situated at a cross-line distance of 180 m and the porosity map at 87 m in depth.The porosity section of Figure12clearly shows high-porosity layers located in the intervals of 85 -87 m and 110 -115 m depth.The production profiles measured in many wells revealed that these layers actually correspond to the major water feeding levels of wellbores.In order to further analyze the spatial distribution of porous bodies and of presumably-conductive flow paths, different cut-off values were applied to the 3D seismic porosity block.Figure13shows several 3D seismic pseudo-porosity blocks, associated with porosity cut-off values of 10 and 30%.The blocks were extracted from the central part of the area.They have a length of 120 m in the in-line direction and a length of 300 m in the cross-line direction.The extracted 3D reservoir volume with a porosity smaller than 10% (Figure13, top right) actually represents the largest fraction of that aquifer, i.e. tight carbonates with a low permeability (less than a millidarcy).This is consistent with the observation of very sparse and channeled flow paths within that aquifer.We may indeed assume that, within that extracted low-porosity volume, the density of conductive (karstic) bodies is too low to ensure a hydraulic communication between wells, because the velocity-to-porosity converted block is derived from a high-resolution 3D

Fig. 14 .
Fig. 14.Diffusivity map from slug tests and seismic porosity map at a depth of 87 m.

Fig. 15 .
Fig. 15.Example of the fitting of fictitious scenario adding sequentially drawdown responses from different pumping tests (Left).Example of hydraulic conductivity field sought by inversion of a fictitious pumping scenario (Right).