FR 2909185 A1 20080530 FR 0610493 A 20061127 10 La présente invention concerne le domaine de l'interprétation stratigraphique de données sismiques. Plus particulièrement, l'invention concerne le domaine des procédés de détection automatique de réflecteurs sismiques montrant des caractéristiques pertinentes pour l'interprétation en stratigraphie sismique, telles que des caractéristiques de convergence de 15 couches sédimentaires. L'imagerie sismique est une méthode d'observation de l'architecture du sous-sol, notamment utilisée dans les domaines de l'exploration pétrolière et de la caractérisation des réservoirs naturels souterrains d'hydrocarbure. Cette technique repose sar l'émission de 20 signaux acoustiques dans le sous-sol et l'enregistrement des signaux propagés et réfléchis sur des réflecteurs acoustiques particuliers. Ces signaux sont traités pour former une image sismique en deux ou trois dimensions. Cette image sismique est constituée d'une série d'enregistrements verticaux appelés "traces sismiques". Ces traces sismiques représentent l'amplitude du signal reçu en fonction du temps. Les signaux enregistrés corrèlent en général 25 d'une trace à l'autre ce qui se traduit sur une image sismique par des lignes sub-horizontales, plus ou moins rectilignes, épaisses et continues appelées linéations, ou encore réflecteurs. Ces lignes représentent les interfaces entre des couches sédimentaires. Une image sismique bidimensionnelle correspond à une section verticale du sous-sol. Un exemple d'image sismique bidimensionnelle est fourni par la figure 2. Les deux dimensions sont la profondeur 30 Z et une direction géographique sub-horizontale X. 1 2909185 A partir de l'observation de l'image sismique, l'interprétateur en stratigraphie sismique dispose de certaines règles, plus ou moins empiriques, pour trouver des indices permettant de reconstituer l'histoire sédimentaire de la zone étudiée, et donc de trouver des zones du sous- 5 sol susceptibles de constituer des pièges pétroliers. Il se base notamment sur l'analyse de la brillance et de la longueur des réflecteurs, sur leur position et leur orientation relativement aux réflecteurs voisins (des réflecteurs parallèles peuvent indiquer un dépôt régulier à grande profondeur d'eau et d'un intérêt pétrolier limité, des réflecteurs convergents organisés en éventail traduisent un basculement progressif du terrain au cours du dépôt et ayant des implications sur la distribution des sédiments ), et sur la manière dont ceux-ci s'interrompent les uns sur les autre. L'analyse de la géométrie des réflecteurs fait par exemple l'objet d'un article fondateur de la stratigraphie sismique: Vail. P.R., R.M. Mitchum, et S. Thompson, III, 1977, Seismic stratigraphy and global changes of sea level, part 4: Seismic stratigraphy and global changes of sea 15 level, in Payton , C.E. (ed.), Seismic Stratigraphy - Applications to Hydrocarbon Exploration: Amer. Assoc. of Petrol. Geologists, Memoir 26, p. 83-97. Cette phase d'interprétation est souvent très délicate et repose entièrement sur l'expertise sédimentologique de l'interprétateur. Il existe très peu d'outils lui permettant de l'assister dans cette tâche, d'autant plus complexe qu'une image sismique comporte souvent plusieurs 20 centaines de réflecteurs en 2D, plusieurs milliers en 3D. État de la technique Plusieurs travaux ont été engagés ces dernières années proposant des méthodes d'analyse d'images sismiques, visant à faciliter la procédure d'interprétation stratigraphique. 25 On peut citer par exemple le document suivant qui fait une synthèse de ces méthodes: S.Copra et K.J.Marfurt, Seismic attribute mapping of structure and stratigraphy , Distinguished Instructor Series, n 9, EAGE, 2006. On peut citer également les travaux du projet TriTex : European research project funded by the European Commiss:on, Directorate-30 General Information Society, "Automated 3D texture content management in largescale data sets", Project No.: IST-1999-20500, Super-final report. 2909185 3 Ces travaux consistent à fournir des méthodes d'analyse d'image, qui transforment une image de départ en une image d'arrivée. En cela, elles ne peuvent constituer qu'un prétraitement avant interprétation. En effet, l'interprétateur sismique doit fournir au final un "pointé" de l'image sismique, c'est-à-dire d'extraire manuellement de celle-ci, par digitalisation, des objets 5 appelés réflecteurs sismiques, exprimés sous forme vectorisée : un réflecteur est vectorisé par une ligne brisée en 2D, ou par une surface discrétisée en 3D. Les méthodes présentées ci-dessus offrent une transformation de l'image de départ qui facilite le pointé ultérieur, mais elles n'effectuent pas ce pointé. Parmi les méthodes d'analyse d'image proposées, certaines sont particulièrement 10 intéressantes pour la stratigraphie sismique, car elles portent sur la détection des convergences dans les images sismiques. On citera par exemple: T.Randen et al, New Seismic Attributes for Automated Stratigraphic Facies Boundary Detection , SEG-98, New Orleans, LA, Expanded Abstracts, September 1998. A.Barnes, Attributes for automating seismic facies analysis , SEG Technical Program, Expanded Abstracts -- 2000 -- pp. 553-556 . Ces méthodes ont en commun la détermination préliminaire d'un "champ de vecteurs directionnels" de l'image à analyser. Il s'agit, en chaque pixel, de déterminer si l'image présente une texture orientée dans un voisinage carré centré sur le pixel considéré, et d'en calculer le(s) vecteur(s) direction. Dans ce dernier calcul, on détermine une unique orientation locale moyenne (Randen, Barnes). Dans les travaux de Barries, la quantification de la convergence des réflecteurs est obtenue en calculant la divergence du champ des vecteurs directionnels. Dans les travaux de Randen, le champ directionnel est soumis à une recherche de lignes de flux dont la densité exprime la divergence ou la convergence des réflecteurs. Ces méthodes travaillent sur des valeurs moyennes calculées dans des voisinages glissants de forme et de taille fixe et prédéterminées, ce qui n'est pas parfaitement adapté à l'interprétation d'images sismiques. En effet, les images sismiques ont la particularité suivante: certains réflecteurs, rectilignes ou non, séparent l'image en deux domaines caractérisés par des textures d'image très différentes dans lesquelles les orientations dominantes de part et d'autre du réflecteur sont différentes. Cela provienr: du fait que les réflecteurs sismiques correspondent à des interfaces sédimentaires représentant des 2909185 4 changements historiques brutaux de conditions sédimentologiques. Par exemple, un réflecteur marque la transition entre un milieu sableux de dunes éoliennes et un milieu marin plus argileux, éventuellement séparés dans le temps par une phase d'érosion. Les méthodes mentionnées ci-dessus ne permettent pas de séparer l'analyse réalisée 5 dans un voisinage placé en fonction de la position d'un réflecteur (complètement au-dessus ou complètement au-dessous d'un réflecteur par exemple), car elles ne connaissent pas la notion de "réflecteur". Elles fournissent une information qui "moyenne" les données provenant d'un voisinage qui peut regrouper des milieux différents et qu'il n'y a pas lieu de caractériser par les mêmes paramètres. 10 L'objet de l'invention est une méthode d'interprétation stratigraphique de données sismiques permettant de caractériser les données par des paramètres propres à chaque réflecteurs sismiques présents dans les données sismiques. L'invention permet notamment de séparer de façon automatique des réflecteurs pour extraire des caractéristiques pertinentes 15 pour l'interprétation stratigraphique, de manière à détecter de façon très précise des agencements particuliers de couches sédimentaires, tels que des agencements de types "onlap" ou "toplap" (convergence de couches sédimentaires liées au type de dépôt sédimentaire). La méthode selon l'invention 20 L'invention concerne une méthode d'interprétation stratigraphique de données sismiques pour déterminer une histoire sédimentaire d'une zone souterraine, dans laquelle lesdites données sont issues de mesures au sein de ladite zone et forment une image multidimensionnelle composée de pixels et comportant des réflecteurs sismiques. La méthode comporte les étapes suivantes : 25 - on digitalise lesdits réflecteurs sismiques à l'aide d'une technique de reconnaissance de forme appliquée à ladite image ; on calcule, pour au moins un réflecteur R, au moins un attribut caractéristique d'une structure de ladite image au voisinage dudit réflecteur, à l'aide d'une technique d'analyse d'image appliquée à au moins une fenêtre d'analyse (FA) dont la forme est dépendante dudit 30 réflecteur R ; 2909185 5 - on sélectionne au moins un réflecteur pertinent au sens de l'interprétation stratigraphique sur la base dudit attribut ; et - on reconstitue l'histoire sédimentaire de ladite zone souterraine à partir d'une interprétation stratigraphique desdits réflecteurs ainsi sélectionnés. Selon l'invention, la technique de reconnaissance de forme peut comporter les étapes suivantes : - on applique une méthode de vectorisation de réflecteurs dans laquelle ladite image est analysée à l'aide d'une fenêtre glissante (FG) centrée successivement sur chaque pixel de 10 l'image, si le pixel placé au centre de la fenêtre d'analyse prend la valeur la plus forte ou la plus faible de la fenêtre d'analyse, celui-ci prend une valeur x, sinon, il prend une valeur différente de x ; - on digitalise lesdits réflecteurs en recherchant des pixels contigus de même valeur et en réalisant une recherche de connectivité ; 15 - on réalise une présélection de réflecteurs de façon à ne garder que les réflecteurs remplissant des critères stratigraphiques définis. On peut, par exemple, réaliser la présélection en ne conservant que les réflecteurs ayant une longueur minimum donnée. Il peut être avantageux d'ajouter une étape de filtrage directionnel. Selon l'invention, pour se conformer à la forme du réflecteur, la fenêtre d'analyse peut 20 comporter N colonnes de hauteur fixe (M1 pixels), ces colonnes étant centrées en hauteur sur le réflecteur R. Cette fenêtre d'analyse peut être translatée verticalement par rapport au réflecteur selon l'une des configurations extrêmes suivantes : totalement en-dessous du réflecteur ou totalement au-dessus du réflecteur. On peut définir une fenêtre d'analyse pour chaque pixel dudit réflecteur R éloigné des 25 extrémités du réflecteur R d'une distance au moins égale à N/2 pixels. Selon l'invention, on peut transformer une fenêtre d'analyse en la déformant par glissement vertical de ses colonnes qui la constituent, de façon à obtenir une fenêtre de forme rectangulaire. La méthode selon l'invention comporte un calcul d'attribut. Ce calcul peut comporter 30 les étapes suivantes : 5 2909185 6 - on calcule au moins un attribut de pixel pour chaque fenêtre d'analyse ; - on trie les pixels de tous les réflecteurs de façon à ne retenir que ceux dont les attributs de pixel respectent des critères choisis par l'interprétateur ; on redéfinit les réflecteurs sur la base des pixels contigus sélectionnés à l'étape 5 précédente; - on calcule des attributs caractérisant un réflecteur à partir des attributs de pixels calculés pour chaque pixel sélectionné ; -on trie les réflecteurs de manière à ne retenir que ceux dont les attributs respectent certains critères choisis par l'interprétateur. Un critère peut consister à défin:Lr une fourchette 10 admissible pour chaque attribut de pixel. Un critère peut également être une appartenance à des classes d'attributs définies à l'aide d'algorithmes de classification statistiques ou par réseaux neuronaux. Enfin, selon l'invention, les attributs de pixels calculés pour chaque fenêtre sont choisis parmi au moins l'un des attributs suivants : 15 - un niveau de gris moyen de ladite fenêtre d'analyse ; - un écart-type des niveaux de gris de ladite fenêtre d'analyse ; un attribut de texture de ladite fenêtre d'analyse ; la direction dominante au sein de ladite fenêtre d'analyse ; - le désordre apparent au sein de ladite fenêtre d'analyse. 20 L'orientation dominante au sein de ladite fenêtre d'analyse, peut être déterminée en réalisant les étapes suivantes : - on construit une première image F, en effectuant pour toute fenêtre d'analyse une transformation fréquentielle complexe suréchantillonnée ; - on construit une seconde image !FI, en calculant le module de ladite première image 25 F ; - on recherche des maxima au sein de cette seconde image, de façon à déterminer deux pics symétriques par rapport au centre de la seconde image IFS définissant ainsi une direction, dont la direction orthogonale correspond à la direction dominante. 2909185 7 D'autres caractéristiques et avantages de la méthode selon l'invention, apparaîtront à la lecture de la description ci-après d'exemples non limitatifs de réalisations, en se référant aux figures annexées et décrites ci-après. 5 Présentation succincte des figures La figure 1 présente l'organigramme illustrant la méthode. La figure 2 présente un exemple de donnée sismique sur laquelle nous avons travaillé. Les figures 3A à 3D représentent schématiquement des étapes du traitement. La figure 4 présente le résultat de la squelettisation. 10 La figure 5 présente la sélection d'un réflecteur particulier. La figure 6 présente un voisinage redressé autour d'un réflecteur particulier. La figure 7 présente le résultat du tri sur l'image complète. Description détaillée de la méthode 15 La méthode selon l'invention permet de réaliser une interprétation stratigraphique de données sismiques, c'est-à-dire de reconstituer l'histoire sédimentaire d'une zone souterraine à partir de données sismiques sous forme d'image 2D ou 3D. Une image sismique est donc constituée d'un ensemble de pixel dont la couleur est fonction de la valeur de la donnée sismique (amplitudes sismiques, impédance,...) en ce pixel 20 auquel correspond une localisation géographique précise au sein du sous-sol. Ces images comportent des lignes sub-horizontales, plus ou moins rectilignes, épaisses et continues appelées réflecteurs sismiques. La méthode permet de caractériser les données sismiques par des paramètres propres à chacun de ces réflecteurs sismiques. L'invention permet notamment de séparer de façon 25 automatique des réflecteurs pour extraire des caractéristiques pertinentes pour l'interprétation stratigraphique. Une caractéristique est pertinente au sens de la stratigraphie sismique lorsqu'elle permet d'apprécier les conditions sédimentaires de dépôt d'une zone souterraine. 2909185 8 Mode de réalisation particulier de l'invention Selon un premier mode de réalisation, l'invention est appliquée à des données sismiques en deux dimensions, comme celles montrées à la figure 2. Les deux dimensions sont la profondeur Z et une direction géographique sub-horizontale X. L'extension à des cubes de 5 données (3D) est explicitée dans un second temps. La méthode comporte une séquence de traitements qui peut se diviser en trois phases : A. Digitalisation des réflecteurs de l'image : phase de pointé automatique ; B. Phase de calculs d'attributs à l'aide de fenêtres adaptées : analyse d'image appliquée à des fenêtres d'analyse dont la forme est contrôlée par les réflecteurs digitalisés 10 automatiquement. C. Phase d'interprétation : tri des réflecteurs sur la base de caractéristiques géométriques diverses, comme par exemple les mesures de convergence calculées à l'étape précédente, de façon à ne retenir que ceux qui sont pertinents au sens de l'interprétation stratigraphique. La succession de ces étapes est détaillée ci-après en se référant à la figure 1, qui illustre la séquence des étapes de la méthode. A) Digitalisation des réflecteurs de l'image : pointé automatique de l'image 20 Filtrage directionnel (étape 1) - FILT On effectue un filtrage directionnel sur les données sismiques afin de :retirer une partie des bruits de l'image que l'on observe, en particulier les bruits directionnels qui peuvent perturber les résultats des étapes suivantes. Le but est donc d'améliorer la continuité apparente des réflecteurs. Cette étape est optionnelle, mais permet d'améliorer la qualité du résultat final. 25 De nombreuses techniques de filtrage sont envisageables pour effectuer cette opération. Un exemple particulièrement avantageux pour la procédure suivant l'invention, est décrit en détail dans les articles suivants : J. Gauthier, L. Duval and J.-C. Pesquet, "A non separable 2D complex modulated lapped transform and its application to seismic data filtering" in Proceedings of Euroropean 30 Signal Processing Conference, September 2005. 15 2909185 et J. Gauthier, L. Duval and J. -C. Pesquet, "Low redundancy oversarnpled lapped transform and application to 3D seismic data filtering" in Proceedings ,9f International Conference on Acoustic Speech and Signal Processing, May 2006 5 Cette méthode se décompose selon les trois étapes suivantes : - on transforme l'image d'origine (donnée sur lesquelles s'applique la. méthode selon l'invention) à l'aide d'une transformée à recouvrement complexe ("lapped transform"). Le principe des transformées à recouvrement est d'extraire des données d'origines des blocs non disjoints, et sur chacun de ces blocs d'appliquer une fenêtre puis une transformée de type 10 fréquentielle complexe redondante. Il s'agit donc d'une transformation locale. - on traite chaque bloc transformé en recherchant la direction dominante (paire de pics dominants) et en retirant tous les coefficients qui ne correspondent pas à des directions proches de celle qui domine. Enfin on effectue un seuillage des coefficients restants pour ne garder que les plus significatifs. 15 -on effectue la transformation inverse qui vérifie les conditions théoriques, dites de "reconstruction parfaite". Cette étape peut-être réalisée en 3D en suivant les même principes que ci dessus. Le détail en est également présenté dans les articles cités précédemment. 20 Recherche et vectorisation des réflecteurs (étape 2) - SQUEL Une méthode de vectorisation des réflecteurs est par exemple proposée dans la demande de brevet EP 1 717 604. Une autre méthode particulièrement bien adaptée pour mettre en oeuvre la procédure selon l'invention est explicitée ci-après. L'image sismique 2D ou 3D est analysée à l'aide d'une fenêtre glissante (FG) se 25 présentant sous la forme d'une colonne de Npixels (de taille par exemple N==5). Cette fenêtre verticale est centrée successivement sur chaque pixel de l'image. Si le pixel analysé, c'est-à-dire si le pixel placé au centre de la fenêtre d'analyse prend la valeur la plus forte de la fenêtre d'analyse, celui-ci prend la valeur 1. Sinon, il prend la valeur 0. Dans l'image obtenue, les pixels de valeurs 1 sont organisés en lignes (en 2D) ou en surfaces (en 3D) ne présentant 30 pas de bifurcations et formant une squelettisation des réflecteurs de l'image sismique. La figure 3A. représente les pixels mis à 1 par un cercle blanc, suivant cette technique. La même 9 2909185 10 technique peut être appliquée en ne retenant que les pixels prenant la valeur ].a plus faible de la fenêtre d'analyse. Les réflecteurs sont ensuite constitués en recherchant les successions de pixels contigus, c'est-à-dire en contact par une face, par un coin ou par une arrête (en 3D). Cela revient donc à 5 extraire des lignes ayant une épaisseur de 1 pixel et formant l'ossature des réflecteurs (Figure 3B). Par leur construction, les lignes (surfaces dans le cas 3D) ainsi extraites présentent la particularité de ne pas être sécantes. Un résultat de cette étape est illustré sur la figure 4. Enfin, les réflecteurs ainsi extraits (digitalisés) sont indexés en faisant une recherche de 8-connectivité en 2D (de 26-connectivité en 3D), c'est-àdire qu'en chaque pixel d'un objet on 10 cherche parmi les 8 (26 en 3D) voisins s'il en existe un qui est non nul, si oui on dit que ces pixels voisins appartiennent au même objet (que l'on appellera dans la suite réflecteur). Présélection des réflecteurs (étape 3) - SLCT La dernière phase du pointé automatique de l'image est constituée par une étape de 15 présélection parmi les réflecteurs digitalisés et indexés lors de l'étape 2, de façon à ne garder que les réflecteurs ayant une longueur significative. Si la longueur des fenêtres d'analyse (définies à l'étape 5) est N, on garde uniquement les réflecteurs dont la longueur est strictement supérieure à N. Ces réflecteurs de longueurs significatives constituent ainsi un ensemble S de réflecteurs à partir desquels les calculs d'attributs sont réalisés. Pour réaliser 20 cette présélection, on peut également choisir tous autres critères stratigraphiques susceptibles d'aider à l'interprétation stratigraphique, autre que la longueur minimale des réflecteurs. B) Calculs d'attributs à l'aide de fenêtres adaptées : analyse d'image La recherche des directions dominantes et les calculs d'attributs dans le voisinage de 25 chaque réflecteur se font de manière itérative, c'est-à-dire que l'on traite les réflecteurs de l'ensemble S un par un. Ce processus est décrit ci-dessous. Initialisation : choix d'un réflecteur de S (étape 4) - CH Après avoir construit un ensemble S de réflecteurs de longueur significative (longueur 30 supérieure à la longueur de la fenêtre d'analyse), la première étape consiste à choisir un 2909185 11 réflecteur R de l'ensemble S qui n'a pas encore été étudié. Puis on applique à ce réflecteur R les étapes 5 à 8. On choisit un autre réflecteur, que l'on traite et l'on réitère cette opération jusqu'à ce que tous les réflecteurs soient traités. 5 Définition de fenêtres d'analyse conformes aux réflecteurs (étape 5) - FEN Le but de cette étape est de sélectionner des pixels en fonction du réflecteur. Pour ce faire on définit pour chaque pixel du réflecteur, des fenêtres d'analyse adaptées au besoin de l'étude stratigraphique. La méthode de construction des réflecteurs assure que chaque réflecteur de S possède au 10 plus un point dans chaque colonne de l'image. Une image est constituée de pixels régulièrement distribués et formant des lignes droites (horizontalement) et des colonnes (verticalement). On peut donc définir des fenêtres d'analyses (FA) qui s'adaptent à la forme du réflecteur R, comme l'illustre la figure 3C. Une fenêtre d'analyse se conforme à un réflecteur lorsque sa forme, son orientation et sa position sont contrôlées par ce réflecteur. Pour se faire, 15 on associe à chaque pixel d'un réflecteur une fenêtre d'analyse. La taille latérale de la fenêtre d'analyse est notée N. Elle dépend de l'échantillonnage temporel des données sismiques. Il est préférable d'englober suffisamment de réflecteurs, ce qui conduit à des valeurs typiques de N de quelques dizaines (16, 32 et 64 par exemple). Pour chaque pixel placé sur le réflecteur R, éloigné des extrémités d'une distance au 20 moins égale à N/2 pixels, on définit une fenêtre d'analyse (FA) centrée en largeur sur ce pixel (PIX). Pour ce faire, on prend en compte les N/2 colonnes de pixels à droite et les N/2 colonnes de pixels à gauche du pixel. La hauteur de chacune de ces fenêtres d'analyse est définie en nombre de pixels M1. Le choix de M1 dépend de la distribution des réflecteurs. Il 25 est préférable que la fenêtre d'analyse recouvre 1 à 3 autres réflecteurs afin de caractériser l'agencement local des réflecteurs entre eux. Ainsi, pour chaque pixel placé sur le réflecteur R, éloigné des extrémités d'une distance au moins égale à N/2 pixels, on définit une fenêtre d'analyse (FA) centrée en largeur sur ce pixel, de taille N colonnes par M1 pixels de haut, se conformant ainsi à la forme du réflecteur, comme l'illustre la figure 3C. 30 Cette fenêtre de taille M1xN peut être centrée en hauteur sur le réflecteur ou décalée verticalement par rapport à celui-ci. L'amplitude du décalage est contrôlée par un paramètre 2909185 12 de décalage Pdec qui correspond à un nombre de pixels. On peut par exemple régler le décalage de manière à placer la fenêtre complètement au-dessus ou au-dessous du réflecteur, afin de caractériser la zone d'image placée au-dessus ou au-dessous du réflecteur. Ces deux solutions extrêmes sont des modes de réalisation particulièrement avantageux de la méthode 5 dans la mesure où ils permettent de caractériser des zones du sous-sol comprises en-dessous ou au-dessus de réflecteurs, et assure une homogénéité des données analysées. Par exemple, la figure 3C illustre la fenêtre d'analyse (FA) du pixel (PIX) placée au centre du réflecteur R constituée de 6 pixels (N=5). Cette fenêtre est entièrement décalée verticalement sous le réflecteur (Pdec 3). 10 Horizontalisation (étape 6) - HOR La fenêtre est ensuite déformée par glissement vertical des colonnes qui la constituent, de manière à obtenir une fenêtre rectangulaire, conformément au schéma de la figure numérotée 3D, ce qui revient à redresser le réflecteur d'intérêt R de telle sorte qu'il soit 15 horizontal. Le résultat de cette étape appliqué à l'exemple de la figure 5 est représenté sur la figure 6. Les fenêtres extraites au-dessus et en-dessous du réflecteur sont notées FA, et FA,. Cette déformation des fenêtres d'analyse permet de les ramener à une forme rectangulaire, adaptée au traitement par des algorithmes d'analyse d'image. 20 Calcul d'attributs caractérisant les fenêtres d'analyse R (étapes 7) - A TTRIB A partir des fenêtres d'analyse, on extrait des données homogènes du point de vue des réflecteurs. C'est-à-dire, que grâce à la méthode selon l'invention, il est possible de sélectionner automatiquement des données le long d'un réflecteur, de part et d'autre, juste au-dessus, ou encore, juste en-dessous de chaque réflecteur d'intérêt (ayant une taille 25 compatible avec la taille de la fenêtre glissante). A partir de ces données propres à chaque réflecteur, on détermine un attribut capable de synthétiser l'information contenue dans ces données, de façon à tirer des conclusions sur le dépôt sédimentaire de la zone souterraine étudiée. Cette étape consiste à déterminer des attributs caractérisant la fenêtre d'analyse horizontalisée ou non. 30 A titre d'exemple, on peut citer différents attributs particulièrement pertinents pour l'interprétation stratigraphique : 2909185 13 - le niveau de gris moyen de l'image constituée par les données sélectionnées par une fenêtre d'analyse FA, comme celles de la figure 6, FAHu ou FAR,. -l'écart-type des niveaux de gris de l'image constituée par les données sélectionnées - des attributs de texture de l'image de la fenêtre d'analyse horizontalisée, notamment son 5 orientation dominante. Ce dernier exemple est décrit en détail ci-après. La direction locale dominante peut être déterminée par différente méthode. La méthode décrite ci-dessous présente l'intérêt d'être adaptée à la nature des données sismiques. 10 Cette analyse se fait localement sur les fenêtres extraites et redressées à l'étape 6 et de taille M1xN. On suit ainsi le réflecteur R redressé pixel après pixel en extrayant à chaque fois une nouvelle fenêtre jusqu'à l'avoir parcouru entièrement. On effectue alors sur chacune de ces fenêtres une transformation fréquentielle complexe suréchantillonnée (par exemple une transformée de type Fourier redondante) ce qui nous conduit a une nouvelle image F. On 15 considère alors l'image formée par le module de l'image transformée, noté classiquement IFS, dans laquelle on va rechercher les maxima. Les données d'entrées étant réelles la transformation complexe conduit à une image IFS symétrique par rapport à son centre. La recherche des maxima conduit donc à deux pics symétriquespar rapport au centre de !FI définissant ainsi une direction. La direction orthogonale à cette direction est appelée direction 20 dominante. Cette méthode de recherche de direction dominante peut être appliquée sur la zone située au-dessus du réflecteur (strictement ou non), sur la zone en-dessous du réflecteur (strictement ou non), ou sur des zones contenant le réflecteur suivant le type de voisinage extrait lors de l'étape précédente. 25 On peut également calculer d'autres attributs permettant d'évaluer le désordre apparent (mélange de directions) dans une image. Par exemple on peut renvoyer la valeur, en module, d'un des pics de la paire symétrique dominante, puis rechercher les valeurs de la deuxième paire de pics symétriques que l'on peut trouver en appliquant un masque directionnel sur la 30 direction principale des pics dans IF]. En réitérant ce procédé, on peut ainsi renvoyer plusieurs valeurs de pics de poids différents qui pourront être utilisées par la suite. 2909185 14 Tous les attributs évoqués ici peuvent être calculés en 3D. La seule différence étant que, dans ce cas, la direction dominante n'est plus représentée par un seul angle comme dans le cas 2D mais pas un couple d'angles définissant un vecteur normal à la surface dominante. 5 Ces attributs, calculés pour chaque fenêtre d'analyse associée à chaque pixel de chaque réflecteur, sont associés aux pixels constituant les réflecteurs. Boucle sur les éléments de S (étapes 8) - RNT Cette partie de la méthode consiste à vérifier si tous les réflecteurs de S ont bien été 10 parcourus. Si c'est le cas (Reps) alors on peut passer à la phase d'interprétation sinon (Rep2) on recommence à partir de l'étape 4. Ainsi, à l'issue de cette étape B) de calculs d'attributs à l'aide de fenêtres d'analyse adaptées, dont la forme est contrôlée par les réflecteurs digitalisés automatiquement, chaque 15 pixel composant les réflecteurs sont informés par des grandeurs pertinentes au sens de l'interprétation en stratigraphie sismique, comme par exemple leur direction relative (convergence) par rapport aux réflecteurs avoisinants. C) Interprétation des attributs : tri des réflecteurs (étape 9) - INTRP 20 Cette phase d'interprétation des attributs, consiste dans un premier temps à déterminer puis affecter des attributs moyens aux réflecteurs, puis à trier les réflecteurs en fonction des valeurs de ces attributs moyens. Ce tri des réflecteurs est réalisé sur la base de caractéristiques géométriques diverses, comme par exemple les mesures de convergence calculées à l'étape précédente, de façon à ne retenir que ceux qui sont pertinents au sens de l'interprétation 25 stratigraphique. Suivant les étapes ci-dessus, on affecte un ensemble de valeurs scalaires ou vectorielles à chaque pixel de chaque réflecteur, que l'on nommera "attributs des pixels". L'étape d'interprétation des attributs consiste à réduire et trier les réflecteurs. Pour ce faire on réalise les étapes suivantes : 2909185 15 - on réduit le nombre d'attributs en les combinant. Un exemple consiste à calculer des continuités ou des désordres apparents à partir des données de valeur des pics que l'on a renvoyées précédemment. En effet en comparant le "poids" de la direction principale et de la direction secondaire on peut en déduire le degré de désordre (mélange de directions) local le 5 long des réflecteurs sélectionnés. Un autre exemple peut être de considérer non plus l'angle au-dessus et en-dessous d'un réflecteur mais la différence des deux. on trie les pixels de tous les réflecteurs et on ne retient que ceux dont les attributs respectent certains critères. Un exemple de critère consiste à définir une fourchette admissible pour chaque attribut de pixel : par exemple, on ne retient que les pixels pour lesquels la 10 direction dominante dans l'image placée au-dessus ou au-dessous dépasse un seuil donné, ce qui traduit une divergence angulaire locale du réflecteur par rapport aux réflecteurs voisins. Un dernier exemple de critère est l'appartenance à des classes d'attributs définies à l'aide d'algorithmes de classification statistiques ou par réseaux neuronaux. - on redéfinit les réflecteurs sur la base des pixels contigus sélectionnés à l'étape 15 précédente. - on calcule des "attributs de réflecteur", par exemple en moyennent les attributs des pixels qui les constituent, ou en analysant la forme du réflecteur : le nombre de pixels composant le réflecteur est un exemple d'attribut de réflecteur qui caractérise sa longueur. La valeur moyenne de niveau de gris est un autre exemple d'attribut, qui caractérise la "brillance" du 20 réflecteur. La valeur moyenne de l'orientation dominante de l'image au-dessus (respectivement au-dessous) du réflecteur est un attribut qui caractérise la divergence angulaire moyenne du réflecteur par rapport aux réflecteurs voisins. - on trie les réflecteurs de manière à ne retenir que ceux dont les attributs respectent certains critères. Un exemple de critère est une fourchette admissible pour chaque attribut de 25 réflecteur. Comme pour le tri des pixels constituant les réflecteurs, des méthodes de tri statistiques ou par réseaux neuronaux peuvent également être utilisées à ce stade. La figure 7 présente le résultat de cette interprétation automatique utilisant les attributs d'angles au-dessus et en-dessous de chaque réflecteur, faisant leur différence, et ne gardant que les angles de plus de 0.2 radian. On voit apparaître sur cette image les réflecteurs 30 correspondants à la zone de convergence. 2909185 16 A partir d'une interprétation stratigraphique de ces réflecteurs minutieusement sélectionnés sur des critères stratigraphiques, il est possible de reconstituer l'histoire sédimentaire de la zone souterraine dont sont extraites les données sismiques. 5 Ainsi, la méthode selon l'invention combine trois approches : la reconnaissance de forme pour rechercher tous les réflecteurs que comporte l'image à analyser, l'analyse d'image pour rechercher à caractériser la texture de l'image au voisinage des réflecteurs détectés dans des fenêtres d'analyse, et enfin, la classification et le tri automatiques de réflecteurs sur la base d' attributs. 10 L'invention fournit ainsi, à partir d'une image sismique quelconque, un ensemble de réflecteurs sélectionnés automatiquement sur des critères stratigraphiques. Cette méthode constitue donc un outil pour l'interprétateur en stratigraphie sismique pour trouver des indices permettant de reconstituer l'histoire sédimentaire de la zone étudiée, et donc de trouver des zones du sous-sol susceptibles de constituer des pièges pétroliers. 15 REVENDICATIONS 1) Méthode d'interprétation stratigraphique de données sismiques pour déterminer une histoire sédimentaire d'une zone souterraine, dans laquelle lesdites données sont issues de mesures au sein de ladite zone et forment une image multidimensionnelle composée de pixels et comportant des réflecteurs sismiques, caractérisée en ce qu'elle comporte les étapes suivantes: - on digitalise lesdits réflecteurs sismiques à l'aide d'une technique de reconnaissance de forme appliquée à ladite image ; - on calcule, pour au moins un réflecteur R, au moins un attribut caractéristique d'une structure de ladite image au voisinage dudit réflecteur, à l'aide d'une technique d'analyse d'image appliquée à au moins une fenêtre d'analyse (FA) dont la forme est dépendante dudit réflecteur R ; - on sélectionne au moins un réflecteur pertinent au sens de l'interprétation stratigraphique sur la base dudit attribut ; et - on reconstitue l'histoire sédimentaire de ladite zone souterraine à partir d'une interprétation stratigraphique desdits réflecteurs ainsi sélectionnés. 2) Méthode d'interprétation stratigraphique selon la revendication 1, dans laquelle ladite technique de reconnaissance de forme comporte les étapes suivantes : on applique une méthode de vectorisation de réflecteurs dans laquelle ladite image est analysée à l'aide d'une fenêtre glissante (FG) centrée successivement sur chaque pixel de l'image, si le pixel placé au centre de la fenêtre d'analyse prend la valeur la plus forte ou la plus faible de la fenêtre d'analyse, celui-ci prend une valeur x, sinon, il prend une valeur différente de x ; - on digitalise lesdits réflecteurs en recherchant des pixels contigus de même valeur et 25 en réalisant une recherche de connectivité ; - on réalise une présélection de réflecteurs de façon à ne garder que les réflecteurs remplissant des critères stratigraphiques définis. 3) Méthode d'interprétation stratigraphique selon la revendication 2, dans laquelle on réalise la présélection en ne conservant que les réflecteurs ayant une longueur minimum 30 donnée. 2909185 18 4) Méthode d'interprétation stratigraphique selon l'une des revendications 2 ou 3, dans laquelle ladite technique de reconnaissance de forme comporte en sus une étape de filtrage directionnel. 5) Méthode d'interprétation stratigraphique selon l'une des revendications 2, 3 ou 4, 5 dans laquelle ladite fenêtre d'analyse comporte N colonnes de hauteur fixe, lesdites colonnes étant centrées en hauteur sur ledit réflecteur R pour se conformer à la forme de ce réflecteur. 6) Méthode d'interprétation stratigraphique selon la revendication 5, dans laquelle la fenêtre d'analyse est translatée verticalement par rapport au réflecteur selon l'une des configurations extrêmes suivantes : totalement en-dessous du réflecteur ou totalement 10 au-dessus du réflecteur. 7) Méthode d'interprétation stratigraphique selon l'une des revendications 5 et 6, dans laquelle on définit une fenêtre d'analyse pour chaque pixel dudit réflecteur R éloigné des extrémités du réflecteur R d'une distance au moins égale à N/2 pixels. 8) Méthode d'interprétation stratigraphique selon l'une des revendications 5, 6 ou 7, 15 dans laquelle la fenêtre d'analyse est déformée par glissement vertical des colonnes qui la constituent, de façon à obtenir une fenêtre de forme rectangulaire. 9) Méthode d'interprétation stratigraphique selon la revendication 1, dans laquelle on calcule ledit attribut en réalisant les étapes suivantes : -on calcule au moins un attribut de pixel pour chaque fenêtre d'analyse ; 20 on trie les pixels de tous les réflecteurs de façon à ne retenir que ceux dont les attributs de pixel respectent des critères choisis par l'interprétateur ; - on redéfinit les réflecteurs sur la base des pixels contigus sélectionnés à l'étape précédente; - on calcule des attributs caractérisant un réflecteur à partir des attributs de pixels calculés pour chaque pixel sélectionné ; - on trie les réflecteurs de manière à ne retenir que ceux dont les attributs respectent certains critères choisis par l'interprétateur. 10) Méthode d'interprétation stratigraphique selon la revendication 9, dans laquelle un critère consiste à définir une fourchette admissible pour chaque attribut de pixel. 2909185 19 11) Méthode d'interprétation stratigraphique selon la revendication 9, dans laquelle un critère est une appartenance à des classes d'attributs définies à l'aide d.'algorithmes de classification statistiques ou par réseaux neuronaux. 12) Méthode d'interprétation stratigraphique selon la revendication 9, dans laquelle 5 lesdits attributs de pixels calculés pour chaque fenêtre sont choisis parmi au moins l'un des attributs suivants : un niveau de gris moyen de ladite fenêtre d'analyse ; un écart-type des niveaux de gris de ladite fenêtre d'analyse ; un attribut de texture de ladite fenêtre d'analyse ; 10 la direction dominante au sein de ladite fenêtre d'analyse ; le désordre apparent au sein de ladite fenêtre d'analyse. 13) Méthode d'interprétation stratigraphique selon la revendication 12, dans laquelle on détermine l'orientation dominante au sein de ladite fenêtre d'analyse, en réalisant les étapes suivantes : 15 - on construit une première image F, en effectuant pour toute fenêtre d'analyse une transformation fréquentielle complexe suréchantillonnée ; - on construit une seconde image IF!, en calculant le module de ladite première image F; - on recherche des maxima au sein de cette seconde image, de façon à déterminer deux 20 pics symétriques par rapport au centre de la seconde image 1FI définissant ainsi une direction, dont la direction orthogonale correspond à la direction dominante.