Monte Carlo Inversion Research Articles

Time series of temperature in Fram Strait determined from the 2008–2009 DAMOCLES acoustic tomography measurements and an ocean model

AbstractA pilot acoustic tomography program in Fram Strait during 2008–2009 measured a year‐long record of acoustic travel times along a 130 km range acoustic path crossing the West Spitsbergen Current. Individual ray arrivals were not observed. Rather, the arrival patterns consisted of a single, stable, broad arrival pulse of about 100 ms duration. Travel time variations of ±0.15 s recorded the vigorous mesoscale environment of the region and the seasonal cycle. To estimate ocean temperature from the tomography data an inverse scheme employed a high‐resolution ocean model for Fram Strait as the reference ocean. The information from the tomographic measurements is primarily average temperature. Estimated temperatures, averaged over 0–1000 m depth and over range, had a mean of 1.11°C and variations of ±0.33°C; the uncertainty of the tomography estimates was about 60m°C. Agreement with an alternate inverse approach based on EOFs and a Markov Chain Monte Carlo inversion scheme relying on a matched‐peak approach was excellent, indicating a robust estimate for ocean temperature. The inverse estimates for average temperature agreed with the equivalent estimates from hydrographic sections obtained along the acoustic path at the start and end of the program. Among other deficiencies, the ocean model greatly underestimated the intensity of the mesoscale fluctuations and exhibited a warm bias of about 0.38°C in section‐averaged temperature. Tomographic measurements in Fram Strait offer unique large‐scale temperature constraints for ocean models through data assimilation. It is anticipated that these constraints will lead to more accurate estimates of the circulation and transports in Fram Strait.

Markov Chain Monte Carlo inversion of temperature and salinity structure of an internal solitary wave packet from marine seismic data

AbstractMarine seismic reflection technique is used to observe the strong ocean dynamic process of nonlinear internal solitary waves (ISWs or solitons) in the near‐surface water. Analysis of ISWs is problematical because of their transient nature and limitations of classical physical oceanography methods. This work explores a Markov Chain Monte Carlo (MCMC) approach to recover the temperature and salinity of ISW field using the seismic reflectivity data and in situ hydrographic data. The MCMC approach is designed to directly sample the posterior probability distributions of temperature and salinity which are the solutions of the system under investigation. The principle improvement is the capability of incorporating uncertainties in observations and prior models which then provide quantified uncertainties in the output model parameters. We tested the MCMC approach on two acoustic reflectivity data sets one synthesized from a CTD cast and the other derived from multichannel seismic reflections. This method finds the solutions faithfully within the significantly narrowed confidence intervals from the provided priors. Combined with a low frequency initial model interpreted from seismic horizons of ISWs, the MCMC method is used to compute the finescale temperature, salinity, acoustic velocity, and density of ISW field. The statistically derived results are equivalent to the conventional linearized inversion method. However, the former provides us the quantified uncertainties of the temperature and salinity along the whole section whilst the latter does not. These results are the first time ISWs have been mapped with sufficient detail for further analysis of their dynamic properties.

Crustal and uppermost mantle structure beneath the United States

AbstractThis paper presents a new model of the shear velocity structure of the crust and uppermost mantle beneath the contiguous U.S. The model is based on more than a decade of USArray Transportable Array (TA) data across the U.S. and derives from a joint Bayesian Monte Carlo inversion of Rayleigh wave group and phase speeds determined from ambient noise and earthquakes, receiver functions, and Rayleigh wave ellipticity (H/V) measurements. Within the Bayesian inverse theoretic framework, a prior distribution of models is posited and a posterior distribution is inferred beneath all of the more than 1800 TA stations across the U.S. The resulting mean and standard deviation of the mean of the posterior distribution at each station summarize the inversion results, which are then interpolated onto a regular 0.25°×0.25° grid across the U.S. to define the final 3‐D model. We present arguments that show that the standard deviation of the posterior distribution overestimates the effect of nonsystematic errors in the final model by a factor of 4–5 and identify uncertainties in density and mantle Q as primary potential sources of remaining systematic error in the final model. The model presents a great many newly resolved structural features across the U.S. that require further analysis and dedicated explication. We highlight here low‐velocity anomalies in the upper mantle that underlie the Appalachians with centers of anomalies in northern Georgia, western Virginia, and, most prominently, New England.

Crustal structure of the Newer Volcanics Province, SE Australia, from ambient noise tomography

Intraplate volcanism is a widespread phenomenon, and is generally regarded to be independent of plate tectonics. The Newer Volcanics Province (NVP) of SE Australia represents the most recent (and arguably still active) intraplate volcanism on the Australian continent, and has been postulated to originate from the combined effect of localized shear flow and edge-driven convection. In this study, we adopt ambient noise tomography and Monte-Carlo inversion methods to construct the first local-scale 3D crustal Vs model of the NVP region with a resolution of ~35km. The model displays distinct crustal velocity features near the eastern and western margins of the NVP, which may point to the existence of a lithosphere-scale plumbing system for the migration of melt associated with the boundary between the Delamerian Orogen and the Lachlan Orogen, which underlies the NVP. In particular, exceptionally high velocities are observed in the middle crust of the Delamerian Orogen, and are best explained by buried magmatic arcs. This interpretation is consistent with a subduction-accretion origin for the Delamerian Orogen. Trans-lithospheric faults might have developed during the accretion processes, providing possible pathways for the migration of NVP magmas to the crust and surface. Our 3D model also images small localized velocity reductions in the lower crust at the region where the two distinct lithospheric units meet. The low velocity zone is spatially correlated with the top of a prominent lithosphere-scale low-resistivity zone (10–30Ωm), which we interpret to represent intruded magmatic sills with small proportions of melt remaining.

Crustal shear velocity structure in the Southern Lau Basin constrained by seafloor compliance

AbstractSeafloor morphology and crustal structure vary significantly in the Lau back‐arc basin, which contains regions of island arc formation, rifting, and seafloor spreading. We analyze seafloor compliance: deformation under long period ocean wave forcing, at 30 ocean bottom seismometers to constrain crustal shear wave velocity structure along and across the Eastern Lau Spreading Center (ELSC). Velocity models obtained through Monte Carlo inversion of compliance data show systematic variation of crustal structure in the basin. Sediment thicknesses range from zero thickness at the ridge axis to 1400 m near the volcanic arc. Sediment thickness increases faster to the east than to the west of the ELSC, suggesting a more abundant source of sediment near the active arc volcanoes. Along the ELSC, upper crustal velocities increase from the south to the north where the ridge has migrated farther away from the volcanic arc front. Along the axial ELSC, compliance analysis did not detect a crustal low‐velocity body, indicating less melt in the ELSC crustal accretion zone compared to the fast spreading East Pacific Rise. Average upper crust shear velocities for the older ELSC crust produced when the ridge was near the volcanic arc are 0.5–0.8 km/s slower than crust produced at the present‐day northern ELSC, consistent with a more porous extrusive layer. Crust in the western Lau Basin, which although thought to have been produced through extension and rifting of old arc crust, is found to have upper crustal velocities similar to older oceanic crust produced at the ELSC.

A Measurement Platform Scheme and Data Post-processing Method for Thermal Conductivity of Li4SiO4 Pebble Bed

The effective thermal conductivity of tritium breeder pebble bed is an important thermal parameter and must be known for the thermo-mechanical design of solid tritium breeder blankets. In order to obtain the parameter, experimental measurement is an effective method. A measurement platform was designed by University of Science and Technology of China for CFETR solid blanket scheme to measure the immediate thermal conductivity data and study the effect of pebble bed temperature, the purge gas pressure and pebble deformation on the thermal conductivity of pebble bed. Measurements were performed based on about 1 mm diameter Li4SiO4 pebbles in the temperature range between 100 and 800 °C, with purge gas pressure ranging from 0.1 to 0.3 MPa. This paper described a measurement platform scheme by thermal probe method. On the other hand, for the sake of increasing the precision of thermal conductivity data transformed from temperature data, some improvements for the data post-processing using Monte Carlo inversion method were made in this paper too.

Image synthesis with graph cuts: a fast model proposal mechanism in probabilistic inversion

Geophysical inversion should ideally produce geologically realistic subsurface models that explain the available data. Multiple-point statistics is a geostatistical approach to construct subsurface models that are consistent with site-specific data, but also display the same type of patterns as those found in a training image. The training image can be seen as a conceptual model of the subsurface and is used as a non-parametric model of spatial variability. Inversion based on multiple-point statistics is challenging due to high nonlinearity and time-consuming geostatistical resimulation steps that are needed to create new model proposals. We propose an entirely new model proposal mechanism for geophysical inversion that is inspired by texture synthesis in computer vision. Instead of resimulating pixels based on higher-order patterns in the training image, we identify a suitable patch of the training image that replace a corresponding patch in the current model without breaking the patterns found in the training image, that is, remaining consistent with the given prior. We consider three cross-hole ground-penetrating radar examples in which the new model proposal mechanism is employed within an extended Metropolis Markov chain Monte Carlo (MCMC) inversion. The model proposal step is about 40 times faster than state-of-the-art multiple-point statistics resimulation techniques, the number of necessary MCMC steps is lower and the quality of the final model realizations is of similar quality. The model proposal mechanism is presently limited to 2-D fields, but the method is general and can be applied to a wide range of subsurface settings and geophysical data types.

The building blocks of Earth and Mars: A close genetic link

The Earth formed in a swarm of Moon- to Mars-sized objects that collided together to build our planet. A large body of work has been dedicated to understanding the Earth's composition as being made of single groups or mixtures of chondrites, however, these models cannot account for the isotopic and elemental characteristics of the Earth. Here, we test mixtures of meteorites, including achondrites, analyzed for seven isotope systems (O, Cr, Ni, Ti, Mo, Ca and Sr), to reproduce the isotope compositions of the Earth and Mars. Our Monte Carlo inversion (a numerical method based on generation of random numbers used to invert multiparameter models) yields a new compositional model where Earth and Mars come almost entirely from the same source material. This finding is in striking agreement with recent planetary formation models in which Earth and Mars formed in a common narrow zone of the protoplanetary disk with Mars being ejected to its current position which prevented further accretion. An important outcome of the model is that a significant mass fraction of the Earth could have been made of volatile depleted and refractory enriched planetary bodies such as angrites (among the oldest known achondrites). This conclusion is also in agreement with new Si isotope data in angrites which suggest that a component of angrites would help explain the difference in δ30Si between the bulk silicate Earth and its building blocks. Our model matches all isotope compositions for both planets, reproduces the volatile element budget of Mars, and accounts for the enrichment in refractory elements of the Earth and Mars compared to chondrites.

Probabilistic 3-D time-lapse inversion of magnetotelluric data: application to an enhanced geothermal system

SUMMARY Surface-based monitoring of mass transfer caused by injections and extractions in deep boreholes is crucial to maximize oil, gas and geothermal production. Inductive electromagnetic methods, such as magnetotellurics, are appealing for these applications due to their large penetration depths and sensitivity to changes in fluid conductivity and fracture connectivity. In this work, we propose a 3-D Markov chain Monte Carlo inversion of time-lapse magnetotelluric data to image mass transfer following a saline fluid injection. The inversion estimates the posterior probability density function of the resulting plume, and thereby quantifies model uncertainty. To decrease computation times, we base the parametrization on a reduced Legendre moment decomposition of the plume. A synthetic test shows that our methodology is effective when the electrical resistivity structure prior to the injection is well known. The centre of mass and spread of the plume are well retrieved. We then apply our inversion strategy to an injection experiment in an enhanced geothermal system at Paralana, South Australia, and compare it to a 3-D deterministic time-lapse inversion. The latter retrieves resistivity changes that are more shallow than the actual injection interval, whereas the probabilistic inversion retrieves plumes that are located at the correct depths and oriented in a preferential north–south direction. To explain the time-lapse data, the inversion requires unrealistically large resistivity changes with respect to the base model. We suggest that this is partly explained by unaccounted subsurface heterogeneities in the base model from which time-lapse changes are inferred.

Linearized two-hydrophone localization of a pulsed acoustic source in the presence of refraction: Theory and simulations

This paper develops an efficient three-dimensional localization method for transient acoustic sources, with uncertainty estimation, based on time differences between direct and surface-reflected arrivals at two hydrophones. The localization method accounts for refraction caused by a depth-dependent sound-speed profile using a ray-theoretic approach for calculating eigenray travel times and partial derivatives. Further, the method provides localization error estimates accounting for uncertainties of the arrival times and hydrophone locations, as well as for depth-dependent uncertainties in the sound-speed profile. In the first of two steps, source depth and range to each hydrophone are estimated using an iterative, linearized Gauss-Markov inversion scheme. In the second step, the estimated source ranges are combined with the hydrophone locations to obtain the source location in the horizontal. Localization performance is analyzed in a simulation study, and the linearized localization estimates and uncertainties are validated by comparison with a fully nonlinear (but numerically intensive) Markov-chain Monte Carlo inversion.

Thin-sheet electromagnetic modeling of magnetovariational data for a regional-scale study

Naturally existing electromagnetic (EM) fields recorded at the surface of the Earth can be used to infer the electrical conductivity distribution of the subsurface. In the magnetovariational sounding (MVS) technique, the transient variations of orthogonal components of the Earth’s magnetic field are measured. In the frequency domain, the magnetic transfer function relates the vertical component to the horizontal components of the magnetic field. This paper describes the thin-sheet modeling of MVS data on a regional scale. The integrated conductivity variations in the thin-sheet model were estimated by applying the Markov chain Monte Carlo (MCMC) inversion algorithm. The application of this method to MVS data from the Finland part of the Fennoscandian Shield has illustrated the utility of both the thin-sheet approximation and MCMC inversion modeling. The conductivity anomalies obtained from this study confirmed the regional-scale geology of the area.

Markov chain Monte Carlo inversion for the rheology of olivine single crystals

AbstractWe present major modifications to the Markov chain Monte Carlo inversion method of Korenaga and Karato (2008), which was developed to analyze rock deformation data and determine a corresponding flow law and its uncertainty. The uncertainties of state variables, e.g., temperature, pressure, and stress, are now taken into account by data randomization, to avoid parameter bias that could be introduced by the original implementation of the cost function. Also, it is now possible to handle a flow law composed of both parallel and sequential deformation mechanisms, by using conjugate gradient search to determine scaling constants. We test the new inversion algorithm extensively using synthetic data as well as the high‐quality experimental data of Bai et al. (1991) on the deformation of olivine single crystals. Our reanalysis of this experimental data set reveals that a commonly adopted value for the stress exponent (∼3.5) is considerably less certain than previously reported, and we offer a detailed account for the validity of our new estimates. The significance of fully reporting parameter uncertainties including covariance is also discussed with a worked example on flow law prediction under geological conditions.

Bayesian linearized two-hydrophone localization of a pulsed acoustic source

A three-dimensional localization method for transient acoustic sources is developed, based on time differences between direct and surface-reflected arrivals at two hydrophones. The method accounts for refraction caused by a depth-dependent sound-speed profile using a ray-theoretic approach, and, further, it provides localization error estimates accounting for uncertainties of the arrival times and hydrophone locations, as well as for depth-dependent uncertainties in the sound-speed profile. In the first of two steps, source depth and range to each hydrophone are estimated using an iterative, linearized Gauss-Markov inversion scheme. In the second step, the estimated source ranges are combined with the hydrophone locations to obtain the source location in the horizontal. Localization performance is analyzed in a simulation study, and the linearized localization estimates and uncertainties are validated by comparison with a fully-nonlinear (but numerically intensive) Markov-chain Monte Carlo inversion. [Wor...

2D characterization of near‐surface : surface‐wave dispersion inversion versus refraction tomography

ABSTRACTThe joint study of pressure (P‐) and shear (S‐) wave velocities ( and ), as well as their ratio (), has been used for many years at large scales but remains marginal in near‐surface applications. For these applications, and are generally retrieved with seismic refraction tomography combining P and SH (shear‐horizontal) waves, thus requiring two separate acquisitions. Surface‐wave prospecting methods are proposed here as an alternative to SH‐wave tomography in order to retrieve pseudo‐2D sections from typical P‐wave shot gathers and assess the applicability of combined P‐wave refraction tomography and surface‐wave dispersion analysis to estimate ratio. We carried out a simultaneous P‐ and surface‐wave survey on a well‐characterized granite‐micaschists contact at Plœmeur hydrological observatory (France), supplemented with an SH‐wave acquisition along the same line in order to compare results obtained from SH‐wave refraction tomography and surface‐wave profiling. Travel‐time tomography was performed with P‐ and SH‐ wave first arrivals observed along the line to retrieve and models. Windowing and stacking techniques were then used to extract evenly spaced dispersion data from P‐wave shot gathers along the line. Successive 1D Monte Carlo inversions of these dispersion data were performed using fixed Vp values extracted from the model and no lateral constraints between two adjacent 1D inversions. The resulting 1D models were then assembled to create a pseudo‐2D section, which appears to be correctly matching the general features observed on the section. If the pseudo‐section is characterized by strong velocity uncertainties in the deepest layers, it provides a more detailed description of the lateral variations in the shallow layers. Theoretical dispersion curves were also computed along the line with both and models. While the dispersion curves computed from models provide results consistent with the coherent maxima observed on dispersion images, dispersion curves computed from models are generally not fitting the observed propagation modes at low frequency. Surface‐wave analysis could therefore improve Vs models both in terms of reliability and ability to describe lateral variations. Finally, we were able to compute VP/VS sections from both and models. The two sections present similar features, but the section obtained from shows a higher lateral resolution and is consistent with the features observed on electrical resistivity tomography, thus validating our approach for retrieving VP/VS ratio from combined P‐wave tomography and surface‐wave profiling.

Design and research on the measurement platform of the effective thermal conductivity for Li4SiO4 and Li2TiO3 pebble bed

China is carrying out the conceptual design of Chinese Fusion Engineering Testing Reactor (CFETR), and the Helium Cooled Pebble Bed (HCPB) blanket concept is one of the main choices for tritium production. Li4SiO4 and Li2TiO3 are the candidate breeder materials for the HCPB blanket concept. In the HCPB blanket, breeding pebbles with the diameter range of 0.6–1.2mm are placed between two plates and the bed shall be cooled. Accordingly, effective thermal conductivity of pebble beds needs to be determined for the heat transfer calculation. Measurements of the heat transfer parameters of Li4SiO4 and Li2TiO3 pebble beds are being performed at the University of Science and Technology of China (USTC). Two measurement methods are being used. One is the steady state method with the use of thermocouples to measure the temperature distribution of the pebble bed. Another is transient thermal probe method using the temperature variation of the thermal probe and Monte Carlo inversion method to calculate the heat transfer parameters of the pebble bed. This paper will report on the progress of these measurement platforms.

Abrupt changes in the rate of Andean Plateau uplift from reversible jump Markov Chain Monte Carlo inversion of river profiles

Knowledge of the rate of surface uplift of the central Andean Plateau provides important boundary conditions for regional geodynamic models and paleoclimate reconstructions. Here we present a fully nonlinear inverse method to extract the rock uplift history of the central Andean Plateau from the Cotahuasi River and its tributaries. Our approach is based on an analytical solution to the stream power model that relates elevation to rock uplift history under the assumptions of constant climatic conditions in space and time and spatially invariant rock uplift. We use a Bayesian framework that allows us to quantify the full state of knowledge of model parameters (i.e., uncertainties, trade-offs). The erosional efficiency defines the landscape response time, and this must be determined using independent data. Therefore, using thermochronometric data from the Cotahuasi–Ocoña Canyon, which record rapid fluvial incision at ~13–10Ma, we calibrate our results to infer a rock uplift history for the north-central Andean Plateau. We infer a pulse of rock uplift between 25 and 10Ma at rates as high as 0.25km/My. However, the rock uplift rate inferred during this pulse is not as high as is predicted if the central Andean Plateau grew as a result of loss of unstable lithospheric mantle. Therefore, our results are more consistent with models requiring gradual uplift of the central Andean Plateau.

Attenuation and velocity structure from diffuse coda waves: Constraints from underground array data

An analysis of coda waves excited in the 0.2–20Hz frequency band and recorded by the underground array Underseis (central Italy) has been performed to constrain both seismic attenuation at regional scale and velocity structure in the Mount Gran Sasso area. Attenuation was estimated with the MLTWA method, and shows a predominance of scattering phenomena over intrinsic absorption. The values of Qi and Qs are compatible with other estimates obtained in similar tectonic environments. Array methods allowed for a detailed study of the propagation characteristics, demonstrating that earthquake coda at frequencies greater than about 6Hz is composed of only body waves. Coherence and spectral characteristics of seismic waves measured along the coda of local and regional earthquakes indicate that the wavefield becomes fully diffuse only in the late coda. The frequency-dependent energy partitioning between horizontal and vertical components has been also estimated and compared with synthetic values computed in a layered half-space under the diffuse field assumption. This comparison confirms that, for frequencies higher than 6Hz, the coda appears as a sum of body waves coming from all directions while, in the low frequency range (0.2–2Hz), the observations can be well explained by a coda wavefield composed of an equipartition mixture of surface and body waves traveling in a multiple-layered medium. A Monte-Carlo inversion has been performed to obtain a set of acceptable velocity models of the upper crust. The present results show that a broadband coda wavefield recorded in an underground environment is useful to constrain both the regional attenuation and the velocity structure of the target area, thereby complementing the results of classical array analysis of the wavefield.

Mapping thunder sources by inverting acoustic and electromagnetic observations

AbstractWe present a new method of locating current flow in lightning strikes by inversion of thunder recordings constrained by Lightning Mapping Array observations. First, radio frequency (RF) pulses are connected to reconstruct conductive channels created by leaders. Then, acoustic signals that would be produced by current flow through each channel are forward modeled. The recorded thunder is considered to consist of a weighted superposition of these acoustic signals. We calculate the posterior distribution of acoustic source energy for each channel with a Markov Chain Monte Carlo inversion that fits power envelopes of modeled and recorded thunder; these results show which parts of the flash carry current and produce thunder. We examine the effects of RF pulse location imprecision and atmospheric winds on quality of results and apply this method to several lightning flashes over the Magdalena Mountains in New Mexico, USA. This method will enable more detailed study of lightning phenomena by allowing researchers to map current flow in addition to leader propagation.

Uncertainty quantification of CO2 saturation estimated from electrical resistance tomography data at the Cranfield site

A parametric bootstrap approach is presented for uncertainty quantification (UQ) of CO2 saturation derived from electrical resistance tomography (ERT) data collected at the Cranfield, Mississippi (USA) carbon sequestration site. There are many sources of uncertainty in ERT-derived CO2 saturation, but we focus on how the ERT observation errors propagate to the estimated CO2 saturation in a nonlinear inversion process. Our UQ approach consists of three steps. We first estimated the observational errors from a large number of reciprocal ERT measurements. The second step was to invert the pre-injection baseline data and the resulting resistivity tomograph was used as the prior information for nonlinear inversion of time-lapse data. We assigned a 3% random noise to the baseline model. Finally, we used a parametric bootstrap method to obtain bootstrap CO2 saturation samples by deterministically solving a nonlinear inverse problem many times with resampled data and resampled baseline models. Then the mean and standard deviation of CO2 saturation were calculated from the bootstrap samples. We found that the maximum standard deviation of CO2 saturation was around 6% with a corresponding maximum saturation of 30% for a data set collected 100 days after injection began. There was no apparent spatial correlation between the mean and standard deviation of CO2 saturation but the standard deviation values increased with time as the saturation increased. The uncertainty in CO2 saturation also depends on the ERT reciprocal error threshold used to identify and remove noisy data and inversion constraints such as temporal roughness. Five hundred realizations requiring 3.5h on a single 12-core node were needed for the nonlinear Monte Carlo inversion to arrive at stationary variances while the Markov Chain Monte Carlo (MCMC) stochastic inverse approach may expend days for a global search. This indicates that UQ of 2D or 3D ERT inverse problems can be performed on a laptop or desktop PC. The parametric bootstrap method provides a promising alternative to the MCMC approach for fast uncertainty quantification of the underlying geophysical inverse problems.

Probabilistic electrical resistivity tomography of a CO2 sequestration analog

article Electrical resistivity tomography (ERT) is a well-established method for geophysical characterization and has shown potential for monitoring geologic CO2 sequestration, due to its sensitivity to electrical resistivity contrasts generated by liquid/gas saturation variability. In contrast to deterministic inversion approaches, probabilistic in- version provides the full posterior probability density function of the saturation field and accounts for the uncer- tainties inherent in the petrophysical parameters relating the resistivity to saturation. In this study, the data are from benchtop ERT experiments conducted during gas injection into a quasi-2D brine-saturated sand chamber with a packingthatmimicsa simpleanticlinal geologicalreservoir. The saturation fields are estimated by Markov chain Monte Carlo inversion of the measured data and compared to independent saturation measurements from light transmission through the chamber. Different model parameterizations are evaluated in terms of the recov- ered saturation and petrophysical parameter values. The saturation field is parameterized (1) in Cartesian coor- dinates, (2) by means of its discrete cosine transform coefficients, and (3) by fixed saturation values in structural elementswhoseshapeandlocationisassumedknownorrepresentedbyanarbitrary GaussianBellstructure.Re- sultsshowthattheestimatedsaturation fieldsare inoverallagreementwith saturationsmeasuredbylighttrans- mission, but differ strongly in terms of parameter estimates, parameter uncertainties and computational intensity. Discretization inthe frequency domain(as inthe discrete cosinetransform parameterization) provides more accurate models at a lower computational cost compared to spatially discretized (Cartesian) models. A priori knowledge about the expected geologic structures allows for non-discretized model descriptions with markedly reduced degrees of freedom. Constraining the solutions to the known injected gas volume improved estimates of saturation and parameter values of the petrophysical relationship.