FR 2635197 A1 19900209 FR 8810719 A 19880805 L'invention concerne une méthode pour obtenir un modèle représentatif d'un milieu hétérogène d'après des mesures indirectes obtenues depuis L'extérieur de celui-ci et d'autres données, notamment un ensemble de données ponctuelles mesurées in situ. La méthode selon l'invention peut convenir dans des domaines très divers selon Le type d'émission empLoyé pour obtenir les images. Elle s'applique notamment dans le domaine de la géophysique et en particulier pour le traitement des données sismiques. Pour les besoins de la prospection sismique par réflexion, des ondes sont émises par une source sismique. Elles se propagent dans les terrains à explorer et l'on capte les réflexions de ces ondes sur Les discontinuités du sous-sol au moyen d'un ensemble de capteurs sismiques répartis Le long de la section de sous-sol à étudier. Les signaux enregistrés en fonction du temps sont généralement échantillonnés sous une forme numérisée. Par un déplacement régulier de la source sismique et de l'ensemble de capteurs et une répétition de cycles d'émission-réception, on enregistre un grand nombre de traces sismiques contenant des informations sur le sous-sol. Divers traitements sont systématiquement appliques aux enregistrements par des calculateurs numériques programmés, pour améliorer leur Lisibilité. Il s'agit en particulier de sommations de traces d'enregistrement par application de la méthode dite de couverture multiple pour améliorer le rapport du signal au bruit ou bien encore de déconvoLutions et de migrations comme il est bien connu des spécialistes, pour améliorer la résolution verticale et la résolution horizontale. -2- En particulier Le traitement en amplitudes préservées permet d'éliminer tous les phénomènes qui ne sont pas liés à La réfLectivité des discontinuités en incidence normale. De cette façon, La section traitée est représentative des interfaces géologiques, Les amplitudes des ondes réfléchies étant alors sensiblement proportionnelles aux coefficients de réflexion correspondant aux discontinuités du sous-sol. Apres le traitement en amplitudes préservées, l'image sismique contient donc des informations Lithologiques. Les coupes sismiques améliorées que l'on obtient après ces traitements sont alors disponibles pour l'interprétation. Il s'agit de transformer chaque coupe sismique en temps en une image lithologique du sous-sol et pour cela donc d'estimer avec suffisamment de précision la valeur de l'impedance acoustique des différentes couches du sous-sol, c'est-à-dire le produit de leur densité par la vitesse de propagation des ondes qui s'y propagent, dont les discontinuités sont responsables des réflexions. Le passage à une représentation lithologique du sous-sol traduisant au mieux les coupes sismiques obtenues, est en général très délicat à résoudre. Une méthode classique permet par une déconvolution dite monocanal et par l'établissement de pseudo-logs d'impédance acoustique, de réaliser cette inversion. Elle consiste essentiellement, l'excitation sismique appliquée au sol étant préalablement estimée, à reconstituer à partir des sismogrammes réels la distribution des valeurs de l'impédance acoustique en fonction de la profondeur mesurée en temps de propagation. Une telle méthode est décrite par exemple dans un article intituLé "Inversion of seismograms and pseudo-velocity logs" in Geophysical Prospecting vol. 25, pp 231-250, 1977. On vérifie le bien-fondé des résultats par référence à d'autres données que l'on possède par ailleurs. En général on dispose - 3 - d'informations sur La position des interfaces et sur certaines caractéristiques du sous-sol, qu'elles résultent d'interprétations antérieures ou aient été obtenues autrement en réalisant des mesures de puits (diagraphies, PSV, etc) dans un ou éventuellement plusieurs forages traversant la partie du sous-sol étudiée, avec différents outils. Les résultats obtenus par ce type de méthode ne sont pas toujours satisfaisants. Ils se revèlent souvent peu conciliables avec les résultats des mesures diagraphiques. La divergence des résultats est généralement attribuée aux incertitudes sur la forme exacte de l'excitation transmise au sous-sol, au bruit de fond important perturbant les signaux sismiques utiles et aussi à leur bande de fréquence relativement étroite. De plus les solutions que l'on peut trouver au problème de l'inversion sont instables. Des distributions d'impédance très différentes peuvent expliquer tout aussi bien une même section sismique. La tâche de l'interprétateur qui doit dériver un modèle géologique d'après la distribution des valeurs d'impédance obtenues par ces inversions, en le conciliant avec les différentes données connues par ailleurs, est donc délicate. On vérifie aussi que dans la pratique, le modèle lithologique obtenu peut être remis en cause. Une approche connue pour remédier aux inexactitudes d'un modèle, consiste à le perfectionner par étapes successives. A partir d'un modèle initial choisi a priori d'après des données résultant d'interprétations et de diagraphies, et à ce stade encore largement hypothétique, on calcule des sismogrammes synthétiques qui constituent la réponse du modèle et on les compare aux sismogrammes réels obtenus. Les écarts constatés traduisant l'imperfection du modèle, on modifie celui-ci. Ce cycle de construction de sismogrammes synthétiques et de comparaison est repris et l'on sélectionne de proche en proche les variations appliquées au modèle de base qui conduisent à une réduction des écarts de manière à se rapprocher le plus possible des sismogrammes reels obtenus. - 4 - IL est égaLement connu d'automatiser Le processus ci-dessus en recherchant La distribution d'impédance acoustique teLLe que Les sismogrammes synthétiques s'ajustent au mieux sur Les enregistrements sismiques. De telles approches en utilisation mono-canal, sont décrites par exemple par Bamberger et aL in Geophysics vol 47, pp 757-770, Mai 1982 ou par Oldenburg et aL in Proceeding of the IEEE, special issue on Seismic Inversion, Mars 1986. L'obtention d'une section en impedance qui soit bien compatible avec les informations géologiques, rend nécessaire que l'on impose une corrélation latérale sur les résultats et donc que l'on trouve des solutions convenables dans un plan plutôt qu'une série de solutions indépendantes en considérant séparément l'inversion des différentes traces. Ces traitements sont multicanaux. Selon le mode de paramétrisation du modèle lithologique, on est conduit à des modes de résolution très différents. L'intervalle (0, X) qui correspond à la longueur du profil sismique que l'on cherche à interpréter, ayant été divisé en une suite d'intervalles disjoints entre des points de contrôle, on peut paramétrer le modèle par la géométrie des couches et en supposant que l'impédance est une fonction constante en profondeur à l'intérieur de chaque couche et continue avec une variation linéaire en x à l'intérieur de chacun desdits intervalles. Une telle méthode est décrite par exemple dans le brevet US 4.679.174. Ce mode de paramétrisation est intéressant s'il permet de réduire le nombre de paramètres à traiter c'est-à-dire si l'épaisseur des différentes couches et les intervalles disjoints ne sont pas choisis trop petits. Mais dans ce cas, la paramétrisation ne reflétera pas la réalité car la distribution des valeurs de l'impédance acoustique à l'intérieur d'une couche géologique est en général bien plus complexe. La méthode selon l'invention permet d'obtenir un modèle optimal à au - 5 moins deux dimensions représentant les variations d'au moins un paramètre physique d'un milieu hétérogène tel que Le sous-sol, responsables des réflexions subies par des ondes émises dans le sous-sol, le modèle s'accordant au mieux avec des enregistrements de signaux reçus par un ensemble de capteurs disposés hors du milieu en réponse aux ondes réfléchies issues du milieu, avec une valeur dudit paramètre physique choisie en au moins un point, obtenue à la suite de mesures ou d'évaluations et avec d'autres informations concernant le type d'hétérogéneéités recherché. Elle comporte: - la construction d'un modèle de référence dans lequel les valeurs dudit paramètre physique s'adaptent aux valeurs choisies et auxdites informations concernant les hétérogénéités, - la construction d'enregistrements synthétiques traduisant la réponse du modèle aux ondes émises depuits l'extérieur et la construction du modèle optimal qui permet un ajustement des enregistrements synthétiques sur les enregistrements de signaux reçus. La méthode selon l'invention est caractérisée en ce que: - la construction du modèle optimal est obtenue en minimisant globalement les écarts entre les valeurs de signaux lus sur lesdits enregistrements de signaux reçus et lesdits enregistrements synthétiques d'une part et les écarts dudit paramètre physique par rapport aux valeurs prises par le même paramètre en tout point du modèle de référence d'autre part, et en ce qu'elle comporte - l'application auxdits écarts d'opérateurs choisis pour modéliser les incertitudes respectivement sur les signaux enregistrés et sur le modèle de référence ainsi que pour imposer au modèle optimal une corrélation entre les valeurs du paramètre aux différents points de - 6- celui-ci. Le milieu hétérogène à modéliser peut être par exempLe une zone du sous-sol, l'ensemble de capteurs étant disposé au-dessus de ladite zone et la valeur du paramètre physique mesurée en au moins un point du milieu est obtenue par des moyens de réception de signaux disposés dans au moins un puits traversant cette zone du sous-sol.- Suivant un mode de réalisation, les signaux enregistrés sont des signaux sismiques renvoyés par les différentes discontinuités du sous-sol, ledit paramètre physique est l'impédance acoustique des différentes couches de terrain, ladite valeur du paramètre est obtenue par des mesures dans au moins un puits, et le modèle de référence représente un modèle stratigraphique dont les valeurs s'accordent aux valeurs d'impédance acoustique obtenues dans chaque puits. Suivant un mode préféré de réalisation, les opérateurs appliqués aux écarts sont des matrices de covariance dont les termes sont choisis en fonction d'une évaluation des incertitudes sur le modèle de référence et les enregistrements de signaux reçus. La méthode selon l'invention, ses effets et ses avantages apparaitront mieux à la lecture de la description ci-après d'un mode de réalisation de la méthode appliquée à la modélisation du sous-sol, donné à titre d'exemple non limitatif, en se référant aux dessins annexes o: - la figure 1 montre un exemple de section sismique à 400 traces (entre 1,05 s et 1,3 s) d'une zone de sous-sol obtenue après des traitements classiques (déconvolution, sommation, etc); - la figure 2 montre un exemple de logs d'impédance acoustique obtenus à partir de diagraphies effectuées dans un puits P1 (A) et dans un puits P4 (B1) traversant ladite zone et un pseudo log (B2) obtenu à l'abscisse du puits (P4) par la méthode selon l'invention; 7 - - la figure 3 montre un modèle de référence construit par la technique de modélisation stratigraphique en prenant comme données les valeurs d'impédance obtenues dans un seul puits P1, Les horizons indiqués par les signaux sismiques étant utilisés pour guider Latéralement L'extrapolation de ces données; - la figure 4 montre une section sismique résiduelle associée au modèle de la figure 3, obtenue par soustraction à la coupe sismique de la figure 1, d'une coupe sismique synthétique reconstituée en prenant pour base le même modèle; - la figure 5 montre un exemple de modèle d'impédance qui peut être obtenu par la méthode selon l'invention si l'interprétateur accorde une confiance très faible au modèle de référence; - la figure 6 montre une section sismique résiduelle correspondant au modèle de la figure 5; - la figure 7 montre un modèle d'impédance optimal obtenu par la méthode selon l'invention avec un choix réaliste d'opérateurs de covariance; - la figure 8 montre une section sismique résiduelle correspondant au modèle optimal de la figure 7; et la figure 9 montre une exemple d'organigramme permettant de construire le modèle optimal recherché. La méthode selon l'invention comporte tout d'abord l'établissement d'un modèle de référence bidimensionnel MR(xx) traduisant les variations d'une caractéristique physique du milieu à partir de données ponctuelles obtenues par des mesures in situ et d'informations stratigraphiques (x représente ici la cordonnée latérale et T la profondeur mesurée en temps de partcours aller retour des ondes depuis -- 8 -- La surface). Dans Le cas de la prospection sismique le modèle représente La distribution de l'impédance acoustique du sous-sol sur une certaine longueur X de profil sismique exploré et une profondeur maximale T (définie en temps). Le modèle à construire doit respecter Les mesures ponctuelles disponibles dans La zone explorée. Elles sont obtenues souvent par des diagraphies de densité et de vitesse qui ont été effectuées dans un ou plusieurs puits. On peut donc connaître en fonction du temps de propagation des ondes entre la source d'émission et les récepteurs sismiques, la distribution des valeurs de l'impédance acoustique des terrains traversés par chaque puits. Des corrections éventuelles ayant été introduites pour mettre chaque distribution d'impédance en accord avec les données sismiques et notamment pour situer correctement la profondeur (en temps) des discontinuités importantes. On peut donc imposer que le modèle bidimensionnel à construire MR corresponde à la verticale de chaque forage, à la suite de valeurs d'impédance qui y ont mesurées. Si I(I l) est le nombre de puits forés aux abscisses x.(i = 1,...,I) pour lesquels on dispose de distributions de valeurs d'impédance Zxi(t), obtenues par des mesures, le modèle doit respecter la condition: MR(Xi,) = Zxi(*) (1) La construction d'un modèle de référence peut se faire par la technique classique de modélisation stratigraphique qui tient compte à la fois des données de diagraphie (relation 1), des informations résultant de l'interprétation structurale des coupes sismiques et de l'interprétation stratigraphique. On observe que les formations souterraines présentent moins d'irrégularités le long des isochrones (horizons indiqués par les signaux) que dans des directions orthogonales. Les isochrones sont - 9 - donc des lignes de corrélation du milieu exploré. La géométrie d'une ligne de corrélation peut être obtenue en pointant certains horizons caractéristiques sur les sections sismiques ayant subi au préalable un traitement classique en géophysique, de déconvolution stratigraphique. L'interprétation structurale permet de diviser la zone du sous-sol explorée en plusieurs parties. Dans chacune d'elles l'interprétateur peut définir les modes de dépôts, à partir de ses connaissances géologiques: débordement (onlap) régression (offlap) ou concordance. Cette définition permet de connaître l'orientation des lignes de corrélation en tout point du milieu. La technique de modélisation stratigraphique consiste à définir des variations linéaires de l'impédance le long de lignes de corrélation d'équation 4 '(x) entre les valeurs connues de Zxi.(4 (xi)) aux X1 "intersections" avec chaque puits. Si l'on ne dispose que de données relatives à un seul puits, on effectue une extrapolation en considérant que le modèle est constant le long des lignes de corrélation. Dans le cas o l'on dispose de données relatives à deux puits par exemple, on peut effectuer entre eux une interpolation linéaire. Le modele de référence de la figure 3 a été construit à partir des données de diagraphie dans un seul puits (P1) en considérant que les paramètres du modèle étaient invariables le long des lignes de correlation. Dans l'exemple étudié, on peut définir un opérateur de modélisation F(Z) permettant de calculer la section sismique synthétique Ssynt(x,t). Cette section représente la réponse sismique de la section en impédance Z(x,'), de la zone du sous-sol explorée sous l'effet des excitations sismiques W (t) supposées connues, appliquées en surface au point d'abscisse x. L'opérateur F a pour expression au point d'abscisse x. L'opérateur F(Z) a pour expression: - 10 - F(Z) = Ssynt(x,t) =5R(xi).Wx(t-1)d c (2) avec d R(x,) = dt log Z (x,i) (3) Les excitations sismiques ou ondelettes W (t) utilisées sont x déterminées préalablement par des opérations classiques de déconvolution stratigraphique, à partir des mesures diagraphiques dans le puits P1 et en supposant ces ondelettes constantes suivant l'axe x. En général le modèle de référence obtenu n'est pas satisfaisant, vis-àvis de l'information sismique, comme le montre la figure 3. Une coupe sismique synthétique Ssynt(x,t) a été produite à partir du modèle stratigraphique puis comparée à la coupe sismique à interpréter Senr.(x,t) . La coupe sismique résiduelle (fig. 4) traduisant les différences entre les deux coupes précédentes montre bien que le modèle stratigraphique ne représente qu'une approximation de la réalité et cela d'autant plus que l'on s'éloigne de l'abscisse ou les données diagraphiques ont été mesurées. La méthode selon l'invention permet de construire dans une seconde étape un modèle optimal représentant au mieux d'après les données disponibles, la distribution M(xI) diun paramètre physique et ici en l'occurrence de l'impédance acoustique. L'interprétateur accordant généralement une grande confiance aux mesures in situ (telles que les données de diagraphie) utilisées pour la construction du modèle de référence, les incertitudes sur celui-ci sont liées à des erreurs sur la géométrie du modèle de référence dans la définition des lignes de corrélation par exemple et à la présence d'hétérogénéités du milieu réel le long de ces lignes. Ces erreurs sont fortement corréLées dans la direction des lignes de corrélation et non corrélées dans les directions orthogonales. On modélise les erreurs sur le modèle de référence par une variation aléatoire gaussienne dont l'opérateur de covariance est choisi à variation exponentielle le long des lignes de corrélation, et diagonal dans les directions orthogonales: Cm(x, x',4,) = V exp.( - covariance est conditionnée par la disponibilité de mesures de puits. La variance V(x,o) et la longueur de corrélation l(x,) peuvent varier latéralement et temporellement. Ces deux paramètres sont donc définis sur l'ensemble du domaine correspondant au modèle de référence. Ainsi la longueur de corrélation doit être nulle ou presque nulle au voisinage d'une discordance. Sinon il se produit un effet de lissage du modele optimal à l'emplacement correspondant. L'interprétateur ayant choisi partout des valeurs convenables pour V(x, covariance Cm (une matrice) à associer au modèle de référence MR(x,). Les mesures indirectes des paramètres physiques de la zone explorée telles que celles obtenues par des opérations de prospection sismique, sont généralement entachées d'erreurs. On traduit ces erreurs par une variable aléatoire gaussienne de moyenne nulle pouvant être corrélée en temps par exemple, le spectre d'amplitude du bruit sismique étant le même que celui du signal, mais non corrélé latéralement (suivant x). Il en résulte que l'opérateur de covariance CS (une matrice également), décrivant les incertitudes sur les mesures indirectes dont on dispose (les données sismiques par exemple) est donc diagonal relativement à la variable x. Il est complètement défini par la - 12 - variance VsCx,t) choisie en chaque point du domaine exploré. L'interprétateur a l'entière responsabilité du choix des variances qui décrivent les incertitudes sur les mesures indirectes dont il dispose. On détermine alors la distribution d'impédance Mopt(x,t) la plus vraisemblable. C'est celle qui, à l'intérieur des limites fixées en x et t, minimise la fonctionnelle suivante: J) = JM) + J s() (5) o: JM(* ) = (M(x,) - Mpr(x,t))T.Cm.(M(x,4) - Mpr(x,)) (6) et (Syt(Xt-Sen(,) 7 J ()=(S (x,t)- Senr(X,t Cs (S xt)) (xt Cs(7) s Synt enr Synt enr -i -1 o les matrices Cm et Cs 1 représentent respectivement les inverses T des matrices C et C choisies et o les termes sous la forme () m s représentent les transposées de matrices. Dans le cas de variations continues la composante JM s'exprime de la manière suivante: JM(xt X)2 + (m(M-M)gra MMpr).u)2] dx dl (6) La détermination du modèle recherché Mopt(x,t) qui minimise la fonctionnelle J peut être obtenue par différentes méthodes connues d'optimisation (méthode dite du gradient ou du gradient conjugué, méthode dite de quasiNewton etc) en accord avec l'organigramme de la figure 9, par des moyens informatiques. L'interprétateur ayant défini un modèle de référence et les valeurs des différents paramètres V et l modélisant les erreurs, la construction du modèle optimal est conduite de façon automatique. Un modèle Mopt(x,), reconstruit en appliquant de la méthode selon - 13 - l'invention et en utilisant le modèle de référence stratigraphique de la figure 3, est représenté à la figure 7. On a pu vérifier sa validité en calculant les sismogrammes synthétiques qui lui sont associés. La coupe sismique résiduelle (fig. 8) traduisant l'écart entre ces sismogrammes synthétiques et les sismogrammes déduits des enregistrements sur le terrain, ne fait apparaître que des "événements" non corrélés latéralement, il ne reste essentiellement que du bruit. La fonctionnelle globale à minimiser J(Z) est la somme d'une fonction d'écart proprement sismique Js(Z) et d'une autre proprement géologique Jm(Z). La première a pour effet de corriger la géométrie du modèle de référence et d'introduire des variations latérales qui expliquent la section sismique déduite des enregistrements. La fonctionnelle géologique J (Z) permet d'imposer au modèle optimal de satisfaire les connaissances géologiques de l'interprétateur. Par un choix approprié de celle-ci, il peut rendre stable la solution du problème inverse qui est fondamentalement sous déterminée, comme il est bien connu des spécialistes. En particulier elle rend possible l'introduction dans la solution, des basses fréquences présentes dans les données de diagraphie mais non dans les données sismiques. Dans le terme géologique, un choix approprié des variations et des longueurs de corrélation a pour effet de structurer latéralement de façon très efficace le modèle obtenu par la mise en oeuvre de la méthode selon l'invention. Ceci a pour effet d'empêcher que la section sismique synthétique déduite du modèle ne s'ajuste sur les bruits incohérents présents dans la section sismique déduite des enregistrements. Il en résulte une bonne élimination du bruit et donc une amélioration de la résolution verticale du sous-sol exploré. On doit aussi souligner que la fonction d'écart géologique introduit un couplage entre les traces adjacentes du modèle d'impédance. L'inversion de la section sismique s'effectue donc globalement et non trace par trace comme on le réalise habituellement (cf. Oldenburg et ai précédemment cité). - 14 - Les relations (4) et (5) montrent que si l'interprétateur accorde peu de confiance aux données sismiques, en choisissant des vaLeurs importantes pour la variance sismique Vs, la fonction sismique J S(Z) joue dans la fonctionnelle J(Z) un rôLe faible par rapport à celui de la fonction géologique JM(C). Le modèle résultant se rapproche aLors du modèle stratigraphique pris comme modèle de référence. Si l'interprétateur au contraire accorde peu de confiance au modèle stratigraphique en choisissant un écart type V grand et une longueur de corrélation L petite, le modèle résultant sera proche d'un modèle obtenu par une classique déconvolution stratigraphique (traitement mono-canal). C'est le cas pour le modèle montré à la figure 5 o l'on a minimisé essentiellement les écarts par rapport aux données sismiques. On constate que le bruit présent sur la section sismique est transféré dans le modèle obtenu, comme le montre La section sismique résiduelle de la figure 6. On voit donc que le modèle optimal est intermédiaire entre un modèle obtenu par la modélisation stratigraphique et un autre résultant d'un traitement du type mono-canal (déconvolution stratigraphique et pseudo-logs) des seules données sismiques. La cohérence structurale du modèle optimal (fig. 7) est meilleure que celle obtenue en utilisant les seules données sismiques, et les bruits incohérents latéralement (suivant x) n'ont pas été interprétés comme le montre la figure 8. Il a été obtenu par extrapolation à partir des mesures effectuées dans un seul puits (P1) et l'on voit sur la figure 2 que les valeurs d'impédance qu'il donne à L'abscisse du puits P4 (courbe B2) correspondent bien à celles qui y ont été effectivement mesurées (courbe B1). La méthode permet donc de prédire des logs à partir de données sismiques. Les opérateurs de covariance Cm et Cs permettent de décrire les incertitudes sur Les mesures sismiques et sur le modèle de référence. Par l'application de la méthode selon l'invention, on peut trouver un - 15 - compromis entre la vraisemblance du modèle de référence et l'ajustement au données sismiques. Si l'on privilégie les données sismiques, la solution tend vers le résultat d'une déconvolution stratigraphique classique tandis que si l'on n'accorde qu'une confiance limitée aux données sismiques, la solution tend vers le modèle de référence. Ainsi, en faisant varier les incertitudes et donc le degré de connaissance sur les données et sur le modèle de référence, L'application de la méthode conduit automatiquement à une grande variété de modèles différents et l'on peut explorer l'ensemble des solutions vraisemblables. Le mode de réalisation qui a été décrit, comporte l'utilisation de valeurs d'impédances résultant des mesures in situ. Si l'on ne dispose pas de telles mesures mais de simples estimations ou évaluations, la méthode selon l'invention peut néanmoins être mise en oeuvre a partir des valeurs de paramètre dont on dispose. Le mode de réalisation qui a été décrit comporte l'utilisation du temps de parcours vertical aller et retour des ondes émises, pour mesurer la profondeur de la zone étudiée. Ceci n'est nullement indispensable. Si l'on dispose de données ou d'estimations concernant la vitesse de propagation C des ondes dans le milieu en fonction de la profondeur, obtenues par quelque procédé qe ce soit, on peut par un changement de variables remplacer la coordonnée par la coordonnée z mesurant la profondeur en mètres, ces deux coordonnées étant reliées par une relation du type: z =. C(x, e) de La méthode selon l'invention a été décrite dans une application à la modélisation de l'impédance acoustique du sous-sol. Il est bien évident qu'elle peut s'appliquer à tout autre paramètre physique tel que par exemple la vitesse de propagation des ondes ou la densité influant surles données sismiques obtenues et que le modèle obtenu - 16 - peut être représenté aussi bien dans le système de coordonnées (x,) que dans le système (x, z). D'une façon plus générale, on ne sortirait pas du cadre de l'invention en appliquant la méthode selon l'invention à la réalisation de tomographies dans des domaines très différents notamment dans le domaine de l'imagerie médicale o du contrôle non destructif. - 17 - R E V E N D I C A T I 0 N S REVENDICATIONS 1. - Méthode pour obtenir un modèle optimal Mopt(x,) à au moins deux 1. Mét Mopt x àaumisdx dimensions représentant les variations d'au moins un paramètre physique d'un milieu hétérogène tel que le sous-sol, responsables des réflexions subies par des ondes émises dans Le sous-sol, le modèle s'accordant au mieux avec des enregistrements (Senr.) de signaux par un ensemble de capteurs disposés hors du milieu en réponse aux ondes réfléchies issues du milieu, avec une valeur dudit paramètre physique choisie en au moins un point et obtenue à la suite de mesures ou d'évaluations, et avec d'autres informations concernant le type d'hétérogénéîtés recherché, la méthode comportant - la réalisation d'un modèle de référence MR(x,Y) dans lequel les valeurs dudit paramètre physique s'adaptent aux valeurs choisies, - la réalisation d'enregistrements synthétiques (Ssynt.) traduisant la réponse du modèle aux ondes émises depuis l'extérieur et - la construction du modèle optimal qui permet un ajustement des enregistrements synthétiques sur les enregistrements de signaux reçus, caractérisée en ce que: - ladite construction du modèle optimal est obtenue en minimisant globalement les écarts entre les valeurs de signaux lus sur lesdits enregistrements de signaux reçus (Senr.) et lesdits enregistrements synthétiques (Ssynt) d'une part et les écarts dudit paramètre physique par rapport aux valeurs prises par le même parametre en tout point du modèle de référence d'autre part, la méthode comportant en outre - l'application auxdits écarts d'opérateurs choisis pour modéliser les incertitudes respectivement sur les signaux enregistrés et sur le - 18 - modèle de référence ainsi que pour imposer au modèLe optimal une corrélation entre les valeurs du paramètre aux différents points de celuici. 2. - Méthode selon la revendication 1, caractérisée en ce que ledit milieu hétérogène est constitué par une zone du sous-sol, l'ensembLe de capteurs est disposé au-dessus de ladite zone et ladite valeur du paramètre est mesurée par des moyens de réception disposés dans au moins un puits traversant ladite zone. 3. - Méthode selon la revendication 1, caractérisée en ce que les signaux enregistrés sont des signaux sismiques renvoyés par les différentes discontinuités du sous-sol et reçus par un ensemble de capteurs disposés au-dessus d'une zone du sous-sol, ledit paramètre physique est l'impédance acoustique Z(x,t) associée aux différentes couches du terrain, ladite valeur du paramètre est mesurée en au moins un puits traversant ladite zone, et le modèle de référence MR représente un modèle stratigraphique dont les valeurs s'accordent aux valeurs d'impédance acoustique obtenues dans chaque puits. 4. - Méthode selon l'une des revendications 1 à 3, caractérisée en ce que les opérateurs sont des matrices de covariance dont les termes sont choisis en fonction d'un évaluation des incertitudes sur le modèle de référence et sur les enregistrements de signaux reçus. 5. - Méthode selon l'une des revendications 1 à 3, caractérisée en ce qu'elle comporte une étape de traitement pour minimiser la somme d'une première fonction JM liée au modèle de référence et d'une seconde fonction J5 liée aux données enregistrées, telles que: TmI( - MR JM = (M - MR) Cm (M - M) M R R et Js = (S - S)T Cs (S -S) synt enr synt enr o - 19 - (M - MR) est un vecteur dont les termes mesurent les écarts estimés du modèle M par rapport au modèle de référence MR, (M - MR)T est le méme vecteur transposé, Cm-1 est l'inverse d'une matrice de covariance traduisant les incertitudes sur Le modèle de référence, (Ssynt - Sen) est un vecteur dont les termes mesurent tes écarts synt enr mesurés entre Les signaux lus respectivement sur Les enregistrements effectués et les enregistrements synthétiques, (S - S)-T est le même vecteur transpose, et synt enr _-1 Cs1 est l'inverse d'une matrice de covariance traduisant les incertitudes sur les enregistrements des signaux captés hors du milieu. 6. - Méthode selon l'une des revendications précédentes, caractérisee en ce que la valeur du paramètre physique en au moins un point est mesurée à L'intérieur du milieu. 7. - Méthode selon l'une des revendications 1 à 5, caractérisée en ce que la valeur du paramètre physique en au moins un point est déterminée par une interprétation de données sismiques. 8. - Méthode selon l'une des revendications 1 à 5, caractérisee en ce que la valeur du paramètre physique en au moins un point, est déterminée par analyse des vitesses de propagation des ondes dans le milieu. 9. - Méthode selon l'une des revendications 1 à 5, caractérisée en ce que la valeur du paramètre physique en au moins un point, est choisie en fonction de la nature du milieu.