画像データ処理
試料の画像データを処理する方法が開示されている。当該方法は、前記試料の少なくとも一部が重なっている複数の空間領域の第1画像と第2画像とのレジストレーションを行う手順、及び、前記のレジストレーションされた画像からのデータを処理することで、前記第1画像と第2画像から得られる画像に加えられる前記試料についての情報を有する合成画像データを取得する手順を有する。
【発明の詳細な説明】
【技術分野】
【０００１】
撮像手法、特に画像データを処理することで、その画像データが取得される物理的試料の特性の知見を得る方法及びシステムに関する。
【背景技術】
【０００２】
地球科学におけるコア材料の画像解析は、基本的には堆積学者及び記載岩石学の分野であった。光学顕微鏡及び電子顕微鏡(SEM)から得られる2次元(2D)画像データは、多孔率、孔のサイズ、グレインサイズ、フローユニット、透水率、速度、及び圧縮率の推定に利用されてきた。従来の記載岩石学の手法によって、鉱物の相情報及び鉱物の相の起源（たとえば摩滅又は自生）を特定することが可能となった。ナノメートル長さのスケールの探索が可能である。
【０００３】
他方、ミクロンスケールの孔構造の詳細な3次元(3D)画像を生成する能力を有するX線断層診断が近年、十分に確立された2D顕微手法にとって有用かつ不可欠なものとして受け入れられてきた。高品質の「ターンキー」X線断層診断システムが利用可能となったことは近年、これらのシステム利用の急速な増大を後押ししてきている。これらのシステムによって、3次元における有孔性材料についての孔/グレイン情報を得ることが可能となる。残念なことに、従来のマイクロCTイメージングでは鉱物の識別が不十分で、従来のマイクロCTイメージングの空間解像度は約1μmに制限される。
【０００４】
コア材料の孔/グレイン/粘土スケールでの特性の評価及び測定への注目が近年高まっている。このスケールでの特性の理解は、石油及びガス産業における応用にとって重要である。コア試料の解析は、カギとなる油層岩石物理学上のフロー特性及び多相フロー特性を生成するのに用いられる。これらの特性は、石油会社が、油田及びガス田の発見、油田及びガス田を生産状態にすること、並びに油田及びガス田の操業の際に直面する高い財務上のリスクを緩和する上で重要である。新たな油田及びガス田を生産状態にすることのコスト又は既存の油田及びガス田の寿命を延ばすことから得られると思われる利益を比較すると、解析それ自体のコストは低い。コアの解析は相反するデータ並びに解釈及び再現の困難なデータをよく与える恐れがあるという事実にもかかわらず、コアの解析は、貯蔵量の推定及び回復速度の予測をする上で業界標準データとなる。そのような困難には、多相フロー特性をよりよく理解するため、基本レベルで解明される必要がある複雑な界面現象が少なくとも部分的に起因している。
【０００５】
コアの実験室内での従来の多相フロー実験において得られた測定結果は、岩石の孔スケール構造と、及び、流体/流体と流体/岩石の相互作用の界面特性の両方の研究に用いられている。多相フローへ適用される孔レベルでの変位のモデル化規則及び方法の開発に大きな関心が存在する。様々な濡れ性のシナリオの下（水で濡れている、混合物で濡れている、大きい/小さい）での、個々の孔内での流体占有についての規則を含む孔ネットワークモデルが、現実の有孔性材料における多相フロー特性の理解を改善させる試みとして発展している（たとえば非特許文献1及び非特許文献2を参照のこと）。今日、様々な流体相の孔ネットワークモデルによる記述を、現実の濡れ性条件下での流体相の分布に関する実験によって得られた孔レベルの情報に対して直接校正するのは不可能である。従って孔のスケールでのモデルを直接孔のレベルで校正するのは不可能である。
【０００６】
CT画像に関する直接的なシミュレーションが、単相の有孔性材料の特性−たとえば透水率、伝導率、及び水銀注入時の毛細管圧力曲線−を予測するのに利用可能であることが示された（たとえば非特許文献3を参照のこと）。
【０００７】
これまでの研究は、マイクロCT画像化実験に基づいて3Dでの孔スケールでの流体分布を特定することが可能であることを示してきた（たとえば非特許文献4を参照のこと）。しかしこれらの研究は、X線CTのビームラインから試料を取り除くことなく、検査中にコア試料をフラッディングさせる必要があることに深刻に悩まされている。このため、3Dにおいて孔スケールで行われる実験は、単純なフラッディング実験に制限されている。
【０００８】
しかもこれらの実験は、長い取得期間にわたってコア材料の位置をミクロンスケールで完全に維持する必要があることで、制約を受けてきた。よって、平衡に達するのにかなりの回数を含む実験（たとえば、油層条件下での定常状態の相対的な透水性測定、天然の原油及びブラインにおけるコアの劣化、有孔性板のフラッディング等）が、数週間またさらには数ヶ月、試料をCT装置で維持するのに必要となるだろう。
【先行技術文献】
【非特許文献】
【０００９】
【非特許文献１】Morrow & Mason, "Recovery of Oil by spontaneous imbibition", Current Opinion in Colloid & Interface Science, Vol. 6, pp. 321-337(2001)
【非特許文献２】H. Behbahani and M. Blunt, "Analysis of Imbibition in mixed wet rocks using pore scale modelling" , SPE 90132, presented at SPE Annual Technical Conference, Houston, 2004
【非特許文献３】C. H. Arns, F. Bauget, A. Ghous, A. Sakellariou, T. J. Senden, A. P. Sheppard, R. M. Sok, W. V. Pinczewski, J. Kelly, and M. A. Knackstedt, "Digital core laboratory: Petrophysical analysis from 3D imaging of reservoir core fragments ", Petrophysics, 46(4), 260-277, 2005.
【非特許文献４】Seright et al., "Characterizing disproportionate permeability reduction using synchrotron X-ray computed tomography", SPE Reservoir Evaluation and Engineering, October 2002, pp. 355-364
【非特許文献５】C. H. Arns, A.P. Sheppard, R. M. Sok and M. A. Knackstedt, NMR petrophysical predictions on digitized core images Petrophysics, 48 (3), 202-221 (2007)
【非特許文献６】M. Knackstedt, C. Arns, A. Sheppard, T. Senden, R. M. Sok, W. Pinczewski and M. Archie 's exponents in complex lithologies derived from, 3D digital core analysis, Ioannidis, Society of Petrophysicists and Well Log Analysts 48th Annual Logging Symposium, Paper Z, June 4-6, 2007, Austin, USA.
【非特許文献７】Morrow, "Wettability and its effect on oil recovery", SPE 21621, published in Journal of Petroleum Technology, December 1990, p.1476
【非特許文献８】Graue et al., "Alteration of wettability and wettability heterogeneity" Journal of Petroleum Science & Engineering, Vol. 33, pp. 3-17
【非特許文献９】M. Prodanovic et al "3D image-based characterization of fluid displacement in a Berea core" Advances in Water Resources 30, 2007, p214
【非特許文献１０】Dautriet et al "Laboratory determination of stress path dependency of directional permeability of Estaillades limestone", Symposium of the Society of Core Analysts, Abu Dhabi, November 2008, paper SCA2008-26
【発明の概要】
【発明が解決しようとする課題】
【００１０】
本発明の目的は、既知の装置の1つ以上の欠点を実質的に克服する、若しくは少なくとも改良すること、又は、有用な代替解決策を提供することである。
【課題を解決するための手段】
【００１１】
本開示の第1態様によると、試料の画像データを処理する方法が供される。当該方法は、前記試料の少なくとも一部が重なっている複数の空間領域の第1画像と第2画像を登録する手順、及び、前記の登録された画像からのデータを処理することで、前記第1像と第2像から得られる画像に加えられる前記試料についての情報を有する合成画像データを取得する手順を有する。
【００１２】
本開示の第2態様によると、前記第1態様の方法に係る工程を実施することが可能なコンピュータプログラムが供される。
【００１３】
本開示の第3態様によると、前記第2態様のコンピュータプログラムを有するコンピュータによる読み取り可能な媒体を有するコンピュータプログラム製品が供される。
【００１４】
本開示のさらなる態様によると、試料の画像データを処理するエレクトロニクスシステムが供される。当該システムは、前記試料の少なくとも部分的に重なっている複数の空間領域の第1画像と第2画像を取得するように備えられた撮像手段、前記の取得された画像を記憶するように備えられた記憶装置、並びに、前記の取得された画像データを取得及び登録することで、前記試料についての情報を有する合成画像データを得る工程であって、前記情報は前記第1画像と第2画像から得られる像に加えられる、工程、を有する。
【００１５】
本開示の他の態様もまた開示されている。
【図面の簡単な説明】
【００１６】
【図１】本発明の一の実施例による方法の概略的なフローダイアグラムを図示している。当該方法は、2D及び3D画像を登録し、その後両画像からの情報を合成することに関する。
【図２Ａ】2D電子顕微鏡像を図示している。
【図２Ｂ】3D断層画像を図示している。
【図３】図2A及び図2Bから得られた登録された画像のアニメーションの最終フレームを図示している。
【図４Ａ】2D光学顕微鏡像を図示している。
【図４Ｂ】3D断層画像を図示している。
【図５】2D顕微鏡像を図示している。
【図６】登録された3D断層画像を図示している。図6の画像に係る解像度は図5の画像に係る解像度とは異なる。
【図７】2D画像から取得された孔構造の情報を3D画像データに合成することによって得られた合成画像を図示している。
【図８】ワーピング変換（warp transformation）が用いられた規則的な正方形グリッドの走査型電子顕微鏡(SEM)像を図示している。
【図９】ワープ変換ステンシルの概略図を示している。
【図１０】ワープ変換において画像の中心に最も近い5つの最大相関値を図示している。
【図１１Ａ】登録前に逆ワープされなかった画像から得られた像を図示している。
【図１１Ｂ】登録前に逆ワープされた画像から得られた像を図示している。
【図１２】本発明の第2実施例による方法の概略的フローダイアグラムを図示している。当該方法は、3D画像同士の位置合わせを行い、その後両画像からの情報を合成することに関する。
【図１３】流体吸入前の試料を図示している。
【図１４】流体吸入後の試料を図示している。
【図１５】試料内部での流体分布を示す動画からのスナップショットである。前記流体分布は3D画像を位置合わせすることによって得られた。
【図１６Ａ】3D相分布を図示している。
【図１６Ｂ】孔のサイズ分布と流体含有量との相関を示している。
【図１７】本発明の第1実施例に係る計算手法のアルゴリズムを図示している。
【図１８】本発明の第2実施例に係る計算手法のアルゴリズムを図示している。
【図１９】図17と図18に示されたアルゴリズムの手順2のフローダイアグラムの詳細を図示している。
【図２０】3D断層画像内部に登録された2D顕微鏡像を図示している。
【図２１】2つの登録された画像中の部位が等しいことを示している。前記2つの登録された画像のうちの一は2DのSEM像で、前記2つの登録された画像のうちの他は3D断層画像から得られた対応部分の薄片である。
【図２２】2DのSEM像と、該2DのSEM像と等価な、登録された3D断層画像から得られた薄片を図示している。
【図２３】3D断層画像から得られた減衰データと2D顕微鏡像から得られた多孔率との間での依存性を表す関数のプロットである。
【図２４Ａ】解像度の異なる2つの画像と、該2つの画像から得られる合成データに基づく画像を図示している。
【図２４Ｂ】図24Aの画像から得られる画像データに基づく、局所的な孔のサイズ分布（Y軸）とX線減衰から得られたミクロ孔マッピングとの相関を示している。
【図２５】3D断層画像から得られた減衰データに対する、2Dミクロ孔から得られた鉱物相の情報のマッピングを図示している。
【図２６】試料から3D画像を取得する手順、前記試料を撮像装置から取り外して、前記試料の実験を行う手順、前記試料を再撮像する手順、及び、前記実験の前後で得られた2つの3D画像を登録することで、前記試料における構造変化又は他の変化を特定する手順の概念を図示している。
【図２７】試料から3D画像を取得する手順、前記試料を撮像装置から取り外して、前記試料の実験を行う手順、前記試料を再撮像する手順、及び、前記実験の前後で得られた2つの3D画像を登録することで、前記試料における構造変化又は他の変化を特定する手順の概念を図示している。
【図２８】カルボン酸で処理された試料の登録3D画像の2つの対応する薄片を図示している。明らかな差異は、試料材料特性への酸の効果に起因する。
【図２９】試料の一連の状態から決まる一連の曲線を図示している。各曲線は、特定の値の水飽和圧力及び毛管圧力によって特徴付けられる。
【図３０】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３１】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３２Ａ】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３２Ｂ】様々な濡れ性条件下での残留炭化水素の小塊のサイズ変化を定量化している。
【図３３】複数の解像度を有する3D画像を登録する応用例を図示している。
【図３４】高解像度の3D画像を低解像度の3D画像へ埋め込む応用例を図示している。
【図３５】直接比較するために配置された異なる解像度の登録3D画像の対応する薄片を図示している。
【図３６Ａ】記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【図３６Ｂ】記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【発明を実施するための形態】
【００１７】
ここで1つ以上の本発明の実施例について図面を参照しながら説明する。
【００１８】
図1乃至図35は、登録された2D画像と3D画像を合成することで得られたデータを処理するための開示された方法に係る様々な態様、及び当該処理方法によって可能となる実際の用途を示している。
【００１９】
A ハードウエアの実装
図36Aと図36Bは一団をなすように、記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【００２０】
図36Aから分かるように、コンピュータシステム3600は、コンピュータモジュール3601、入力装置−たとえばキーボード3602、マウスポインタ装置3603、スキャナ3626、及びマイクロホン3680−、及び、出力装置−プリンタ3615、ディスプレイ装置3614、及びラウドスピーカー3617を含む−によって形成される。1つ以上の撮像装置3627もまたコンピュータモジュール3601と直接接続して良い。これらの装置は、後述の方法によって処理される画像データを供する。
【００２１】
外部の変調器−復調器（モデム）トランシーバ装置3616は、接続3621を経由してコミュニケーションネットワーク3620とのやり取りを行うため、コンピュータモジュール3601によって使用されて良い。ネットワーク3620は広域ネットワーク(WAN)−たとえばインターネット又はプライベートWAN−であって良い。接続3621が電話線である場合、モデム3616は従来の「ダイアルアップ」モデムであって良い。あるいはその代わりに、接続3621が高容量（たとえばケーブル）接続である場合、モデム3616はブロードバンドモデムであって良い。ネットワーク3620へのワイヤレス接続用にワイヤレスモデムが用いられても良い。
【００２２】
コンピュータモデム3601は一般的に、少なくとも1つのプロセッサ、及び、たとえば半導体ランダムアクセスメモリ(RAM)及び半導体読み出し専用メモリ(ROM)から構成されるメモリユニット3606を有する。図36Aに図示されたコンピュータシステムは2つのプロセッサ3605Aと3605Bを有する。2つのプロセッサ3605Aと3605Bは、後述する方法の手順の少なくとも一部において実装される並列処理を可能にする。モジュール3601はまた多数の入出力(I/O)インターフェースをも有する。前記多数の入出力(I/O)インターフェースは、ビデオディスプレイ3614、ラウドスピーカー3617、及びマイクロホン3680と結合するオーディオ−ビデオインターフェース3607、キーボード3602、マウス3603、スキャナ3626、撮像装置3627、及び任意のジョイスティック（図示されていない）用のI/Oインターフェース3613、並びに、外部モデム3616及びプリンタ3615用のインターフェース3608を有する。一部の実施例では、モデム3616はコンピュータモジュール3601内部−たとえばインターフェース3608内部−に組み込まれて良い。コンピュータモジュール3601はまたローカルネットインターフェース3611をも有する。ローカルネットインターフェース3611は、接続3623を介して、コンピュータシステム3600を、ローカルエリアネットワークとして知られているローカルコンピュータネットワーク3622へ接続することを可能にする。また図示されているように、ローカルネットワーク3622はまた、接続3624を介して広域ネットワーク3620へも接続して良い。広域ネットワーク3620は一般的には、所謂「ファイアーウォール」デバイス又は類似の機能を有するデバイスを含む。インターフェース3611は、Ethernet（商標）回路カード、Bluetooth（商標）ワイヤレス構成、又はIEEE802.11ワイヤレス構成によって形成されて良い。
【００２３】
インターフェース3608と3613はシリアル接続及び/又はパラレル接続が可能である。シリアル接続は一般的に、ユニバーサルシリアルバス(USB)標準規格に従って実装されて良く、かつ対応するUSBコネクタ（図示されていない）を有する。他の記憶装置−たとえばフロッピーディスクドライブ及び磁気テープ装置（図示されていない）−が用いられても良い。光学ディスクドライブ3612は一般的には、不揮発性データ源として機能するように供される。持ち運び可能なメモリ装置−たとえば光学ディスク（たとえばCD-ROM,DVD）、USB-RAM、及びフロッピーディスク−が、システム3600への適切なデータ源として用いられて良い。
【００２４】
コンピュータモジュール3601の構成部品3605乃至3613は一般的に、相互接続バス3604を介して、当業者に既知のコンピュータシステム3600が従来の動作モードとなるようにやり取りを行う。上述の構成を実施することが可能なコンピュータの例には、IBM-PC、サンマイクロシステム(Sun)のスパークステーション(Sparkstations)、アップル(Apple)社のマッキントッシュ(Macintosh)（商標）又はこれらの進化型が含まれる。
【００２５】
上述の画像処理方法はコンピュータシステム3600を用いて実装されて良い。これから説明する図1乃至12の手順及び図17乃至19の計算アルゴリズムは、コンピュータシステム3600内部で実行可能な1つ以上のソフトウエアアプリケーションプログラム3633として実装されて良い。特に、上述の方法の手順は、コンピュータシステム3600内部で実行されるソフトウエア3633内の命令によって行われる。ソフトウエアの命令3631は1つ以上のコードモジュールとして生成されて良い。各コードモジュールは特定のタスクを実行することを意図している。そのソフトウエアはまた2つの別個の部分に分割されて良い。第1部分と対応するコードモジュールは当該方法の手順を実行する。第2部分と対応するコードモジュールは前記第1部分とユーザーとの間のユーザーインターフェースを管理する。
【００２６】
ソフトウエア3633は一般的に、コンピュータによる読み取り可能な媒体からコンピュータシステム3600へロードされ、続いて典型的には、図36Aに図示されているようにHDD3610に記憶されるか、又はメモリ3606に記憶される。その後ソフトウエア3633はコンピュータシステム3600によって実行可能となる。一部の例においては、アプリケーションプログラム3633は、ユーザーへ供給され、1つ以上のCD-ROM3625上でエンコードされ、かつメモリ3610又は3606に記憶される前に、対応するドライブ3612を介して読み出されて良い。あるいはその代わりにソフトウエア3633は、ネットワーク3620又は3622からコンピュータシステム3600によって読み出されて良いし、又は他のコンピュータによる読み取り可能な媒体からコンピュータシステム3600へロードされても良い。コンピュータによる読み取り可能な媒体とは、実行及び/又は処理のためコンピュータシステム3600へ命令及び/又はデータを供することに関与する如何なる記憶媒体をも指称する。そのような記憶媒体の例には、フロッピーディスク、磁気テープ、CD-ROM、ハードディスクドライブ、ROM、集積回路、USBメモリ、磁気光学ディスク、又はPCMCIAカードのようなコンピュータによる読み取り可能な媒体等が含まれる。これらはコンピュータモジュール3601の内部又は外部であるかを問わない。コンピュータモジュール3601へのソフトウエア、アプリケーションプログラム、命令、及び/又はデータの提供にも関与することが可能なコンピュータによる伝送可能な媒体の例には、高周波又は赤外伝送チャネル及び他のコンピュータ若しくはネットワーク装置へのネットワーク接続、並びに、eメール伝送及びウエブサイトに記録された情報を含むインターネット及びイントラネット等が含まれる。
【００２７】
上述のアプリケーションプログラム3633の第2部分及び対応するコードモジュールは、ディスプレイ3614上に描画又は表示される1つ以上のグラフィカルユーザーインターフェース(GUI)を実装するように実行されて良い。典型的にはキーボード3602とマウス3603の操作を介して、コンピュータシステム3600とアプリケーションのユーザーは、機能が一致するように前記（複数の）GUIに係るアプリケーションへの制御コマンド及び/又は入力を供するようにしてインターフェースを制御することができる。機能が一致するような他の形式のユーザーインターフェースが実装されても良い。そのようなインターフェースとはたとえば、ラウドスピーカー3617を介したスピーチプロンプト出力を利用するオーディオインターフェース及びマイクロホン3680を介したユーザー音声コマンド入力である。
【００２８】
図36Bは、プロセッサ3605A又は3605Bのうちのいずれか1つ及び「メモリ」3634の詳細な概略ブロック図である。メモリ3634は、図36Aのコンピュータモジュール3601によってアクセス可能な全メモリ装置の論理的集合を表す（HDD3610及び半導体メモリ3606も含まれる）。
【００２９】
コンピュータモジュール3601は最初電源がオンにされるとき、ハードウエアの初期化と動作チェック(POST)プログラム3650が実行される。POSTプログラム3650は一般的には半導体メモリ3606のROM3649内に記憶される。ハードウエア装置−たとえばROM3649−のような恒久的に記憶されるプログラムはファームウエアと呼ばれることがある。POSTプログラム3650は、適切に機能することを保証するためにコンピュータモジュール3601内部でハードウエアを点検し、かつ典型的には、プロセッサ3605A又は3605Bのうちの対応する1つ、メモリ3609又は3606、及び、校正操作を行うために一般的にはROM3649内に記憶される基本入出力システムソフトウエア(BIOS)モジュール3651をチェックする。一旦POSTプログラム3650が問題なく動作すると、BIOS3651はハードディスクドライブ3610を起動する。ハードディスクドライブ3610が起動することで、該ハードディスクドライブ3610内に存在するブートストラップローダープログラム3652が、プロセッサ3605A又は3605Bのうちの対応する1つによって実行される。これにより、OS3653はRAMメモリ3606へロードされ、RAMメモリ3606上でOS3653は動作を開始する。OS3653は、様々な高次の機能を実現するための、プロセッサ3605A又は3605Bのうちの対応する1つによって実行可能なシステムレベルのアプリケーションである。様々な高次の機能には、プロセッサ管理、メモリ管理、デバイス管理、記憶装置管理、ソフトウエアアプリケーション管理、及び汎用ユーザーインターフェース管理が含まれる。
【００３０】
OS3653は、コンピュータモジュール3601上で動作するプロセス又はアプリケーションの各々が十分なメモリを有することを保証するため、他のプロセスに割り当てられたメモリと相反することなく実行できるようにメモリ3609又は3606を管理する。さらに各プロセスが効率的に動作できるように、システム3600において利用可能な様々な型のメモリが適切に用いられなければならない。従って、メモリの一団3634は、どのようにしてメモリの特別な部分が割り当てられるのか（指定されない場合）を示すことを意図しておらず、むしろコンピュータシステム3600によってアクセス可能なメモリの一般的な視点及びそのようなメモリがどのようにして用いられるのかを示すことを意図している。
【００３１】
プロセッサ3605Aと3605Bの各々は多数の機能モジュールを有する。前記多数の機能モジュールには、制御ユニット3639、算術論理ユニット(ALU)3640、及び局所又は内部メモリ3648−キャッシュメモリと呼ばれることもある−が含まれる。キャッシュメモリ3648は一般的に、レジスタ部内に多数の記憶レジスタ3644-3646を有する。1つ以上の内部バス3641はこれらの機能的モジュールと機能的に相互接続する。プロセッサ3605Aと3605Bの各々は一般的に、接続3618を用いて、システムバス3604を経由して外部装置とのやり取りを行う1つ以上のインターフェース3642をも有する。
【００３２】
アプリケーションプログラム3633は一連の命令3631を有する。一連の命令3631は、条件的分岐及び巡回命令を有して良い。プログラム3633はまた、プログラム3633の実行に用いられるデータをも有して良い。命令3631及びデータ3632は、メモリ位置3628-3630及び3635-3637にそれぞれ記憶される。命令3631とメモリ位置3628-3630の相対的なサイズに依存して、メモリ位置3630内に図示された命令によって示されているように、特定の命令が単一のメモリ位置に記憶されて良い。あるいはその代わりに、命令は多数の部分に分割されて良い。メモリ位置3628-3629内に図示されている命令セグメントによって示されているように、各部分は各独立したメモリ位置に記憶される。
【００３３】
一般的には、プロセッサ3605Aと3605Bには、内部で実行される1組の命令が与えられる。続いてプロセッサ3605Aと3605Bは、他の組の命令の実行によって応答するまで、後続の入力を待つ。各入力は多数の入力源のうちの1つ以上から供されて良い。前記入力には、1つ以上の入力装置3602と3603によって生成されるデータ、ネットワーク3620と3622のうちの1つを経由して外部から受信されるデータ、記憶装置3606と3609のうちの1つから取得されるデータ、又は、対応する読み取り装置3612へ挿入された記憶媒体3625から取得されるデータが含まれる。1組の命令が実行される結果、場合によっては、データが出力されることがある。実行には、メモリへのデータ又は変数の記憶もが含まれて良い。
【００３４】
本明細書において以降に開示される計算機の構成は、メモリ3634内の対応するメモリ位置3655-3658内に記憶される入力変数3654を用いる。実行された計算は、メモリ3634内の対応するメモリ位置3662-3665内に記憶される出力変数3661を生成する。中間変数が、メモリ位置3659、3660、3666、及び3667内に記憶されても良い。
【００３５】
動作するプロセッサ3605A並びに/又はプロセッサ3605Bのレジスタ部3644-3646、算術論理ユニット(ALU)3640、及び制御ユニット3639は、協働することで、プログラム3633を構成する命令の組に含まれる命令毎に「取得、復号、及び実行」のサイクルを行う必要があるミクロ動作のシーケンスを行う。各取得、復号、及び実行サイクルは、
(a) メモリ位置3628からの命令3631の取得又は読み取りを行う取得動作、
(b) 制御ユニットがどの命令が取得されたのかを判断する復号動作、
(c) 制御ユニット3639及び/又はALU3640が前記命令を実行する実行動作、
を有する。
【００３６】
その後、次の命令のためにさらなる取得、復号、及び実行サイクルが実行されて良い。同様に、制御ユニット3639がメモリ位置3632へ値を記憶又は書き込みするための記憶サイクルも実行されて良い。
【００３７】
その後、図17-19のアルゴリズムに含まれる各手順又は下位の処理が、プログラム3633の1つ以上の部分に関連付けられ、かつ、プロセッサ3605A（及び/又は3605B）内のレジスタ部3644-3647、ALU3640、及び制御ユニット3639によって実行される。プロセッサ3605A（及び/又は3605B）内のレジスタ部3644-3647、ALU3640、及び制御ユニット3639は、協働することで、プログラム3633を構成する命令の組に含まれる命令毎に「取得、復号、及び実行」のサイクルを行う必要があるミクロ動作のシーケンスを行う。
【００３８】
あるいはその代わりに、画像処理の開示された方法は、様々な計算ルーチンの機能又は下位の機能を実行する専用ハードウエア−たとえば1つ以上の集積回路−で実装されても良い。そのような専用ハードウエアには、グラフィックプロセッサ、デジタル信号プロセッサ、又は1つ以上のマイクロプロセッサと付属のメモリが含まれて良い。
【００３９】
2つの画像のデータの組が有する情報の比較及び合成を可能にするため、各対応する画像データの組のうちの少なくとも1つは、2つの組が1つの座標系を表すように変換されなければならない。光イメージングの分野では、この変換は「レジストレーション（画像の位置合わせ）」と呼ばれ、かつ、1つの座標系を表すように「位置合わせ」される画像は「レジストレーションされた」画像と呼ばれる。
【００４０】
B レジストレーションされた2D画像と3D画像の処理
以降の開示は、3D撮像−たとえば従来のマイクロCT−の能力を拡張し、かつ2D顕微手法の範囲で得られる相補的データ及び/又は高解像度データを用いることによる、3D画像データの校正を可能にする方法に関する。そのようにして得られた合成データの一の用途は、コア材料/試料の地質学的評価である。具体的には、開示された方法は、3D画像上及び該3D画像内部で2D画像の画像データのレジストレーションに関する。
【００４１】
3D画像体積内部に存在する2D顕微鏡画像のピクセル−ボクセル識別を3D画像と同程度の解像度で行うことで、3D画像中のX線グレースケールレベルを顕微手法から得られたある範囲の特性にマッピングすることが可能となる。また前記ピクセル−ボクセル識別は、鉱物学、元素組成、表面エネルギー、及び音響/弾性特性のうちの少なくとも1つを含む特性の研究に用いられる。従来のX線マイクロCTが約1μmの解像度を有する一方で、顕微手法は、部位及び特性の探索をナノスケールで行うことを可能にする。試料内部の少なくとも一部が重なっている複数の空間領域の3D画像に対して2D画像を位置合わせして、該2D画像を3D画像へマッピングする能力は、現時点においてX線CT解析では微細すぎると考えられている特性についての情報の取得を可能にする。レジストレーションされ、かつ処理された画像は、任意に配向し、スケールを変換し、かつ相互に並進して良い。さらに2つの3D画像の処理のレジストレーションが行われる場合−この場合2つのレジストレーションされた3D画像が任意に配向し、スケールを変換し、かつ相互に並進して良い−にも、開示されたこの有利な能力はまた適用可能である。2D画像と3D画像との間でのレジストレーションの場合では、部分的な重なりは、3Dイメージング手法によってイメージングされる3D体積と、2Dイメージングによって試料上又は該試料内部にイメージングされる任意に配向した面との間での線又は面の重なりであることに留意して欲しい。3D-3Dレジストレーションの場合では、部分的な重なりは一般的には、1つ以上の各対応するイメージング手法によってイメージングされる、試料の2つ以上の3D領域間での空間的な重なりである。
【００４２】
開示された方法によって実行されるマッピングは、ある範囲の2D顕微手法によって得られた従来の堆積データ（たとえばコアの記載、記載岩石学、地質学上の外見、地質学上の岩石の種類、鉱物の相、鉱物の起源、弾性/音響特性、表面化学に関連するデータ）を、マイクロCT解析によって得られる3D画像データ（たとえば減衰、グレイン形状、グレインの接続性、孔の接続性、孔の形状に関連するデータ）に相関させることを可能にする。顕微手法では、ナノメートルでの探索が可能である。これはマイクロCTの解像度よりも小さい。続いて得られた高解像度データは、3DマイクロCT画像からのX線減衰データに直接相関付けられて良い。よって従来の堆積学的情報は、油層岩石物理学的特性（透水性、排水毛管圧力、音響及び地震特性、抵抗力、NMR）及び油層の設計に係る特性（たとえば相対的透水性、残留飽和度、吸水毛管圧力）の予測に関連するデータに相関付けることができる。
【００４３】
3D画像に対する2D画像のレジストレーションから得られる合成情報
当該方法の詳細が図1に図示されている。一旦有孔性材料が取り出されると、手順102では、マイクロ断層撮像によって3D画像が得られる。試料から得られた1つ以上の2D部分又は露出した平面は、手順104において、顕微鏡観察/分光用に収集及び調製される。コア材料の残りは品質制御実験（たとえばヘリウム及び水銀注入による多孔度測定）に用いられて良い。
【００４４】
顕微鏡観察/分光は、手順106において、以下の手法のうちの1つ以上によって行われる。
− 走査電子顕微鏡(SEM)：2次電子及び後方散乱電子
− 光学顕微鏡
− 走査音波顕微鏡
− レーザー共焦点顕微鏡
− 集束イオンビーム−走査電子顕微鏡(FIBSEM)
− エネルギー分散分光(EDAX)
− 2次イオン質量分光(SIMS)
− 分光学的手法：赤外、UV−可視、及びラマン分光が含まれるが、これらに限定されるわけではない
− 他の顕微又は分光手法
得られた顕微/分光2D画像データは、従来手法で処理されることで、ある範囲の長さスケール（nm〜μm）にわたって、2Dにおいて試料の様々な構造の特性が取り出される（手順110）。この範囲であれば、たとえば鉱物学、表面特性、局所的な表面相互作用、音響特性等の特性の解析は可能である。例として、光学顕微鏡による標準的な記載岩石学的解析は、以下のうちの少なくとも1つに係る特性の測定を可能にする。
− 有孔性
− 鉱物（たとえば石英、長石、方解石、粘土、黄鉄鉱等）
− グレイン境界、セメント相の特定
− 地質学的外見の種類の特定
− 地質学上の岩石の種類
当該方法の手順108は、3D画像上の1つ以上の2D画像のピクセル−ボクセルレジストレーションを実行する。手順110は従来の2D画像解析（たとえば分光写真解析）を実行する。これにより、当該方法の手順112において、定量的な合成データを得ることが可能となる。手順112は、1種類の手法によって得られたデータに値を追加することが可能である。合成された情報の一の用途は、高解像度の高品質2D顕微鏡解析データとの比較によって、3D画像データ上での直接的な品質制御を行うことである。合成データはまた、高解像度2D顕微鏡画像から得られた孔及び鉱物相の情報を、サブCT解像度又は従来のCTの解像度を有する3D画像から得られたデータに直接相関付けるのに用いることもできる。そのような相関は以下のうちの1つに用いられて良い。
− マイクロ多孔性測定の直接検査
− サブCT解像度スケールから得られたグレイン/孔のサイズ情報のマイクロ多孔性データへの追加
− 2D顕微手法及び分光手法によって得られた情報（たとえば鉱物学、膠着(cementation)、続成作用の履歴、高解像度（マイクロ）多孔性、実効透水率、グレインサイズ、サブミクロン解像度での孔のサイズ、音響特性、及び表面化学特性）の、マイクロCTから得られた3D画像データへの追加
また上述のデータ補正を用いることによって、ある範囲の物理的特性を任意の手順114においてシミュレーションすることができる。そのような物理的特性には、様々な流体相の相対透水率、誘電応答、NMR応答、及び音響地震応答が含まれる。典型的な適用可能なシミュレーションについての情報は、特許文献3、特許文献5、及び特許文献6で見つけることができる。
【００４５】
当該手法の応用例
1. 一例では、3D画像データの直接的な品質制御は、該3D画像データを標準的な2D顕微鏡解析から得られたデータと比較することによって実現されて良い。図2及び図3は、岩石試料のピクセル−ボクセルレジストレーションを表している。前記岩石試料の孔のほとんどはマイクロCTデータの解像度で十分識別されうる。図2AはSEM画像を図示している。図2Bは、手順108（図1）に示されているように、開示された方法によってレジストレーションがなされた各対応する断層薄片を図示している。2D画像と3D画像との直接的な対応が存在する。如何なるばらつきも顕微鏡解析の試料調製に起因する。
【００４６】
2. 他の応用では、3D画像データからの減衰情報と、2D顕微鏡解析からの改善された情報とが直接的に相関付けられてよい。図4Aは、砂岩石試料の2D光学顕微鏡画像の一例を図示している。鉱物学上の広範囲にわたるばらつきを観測することができる。鉱物相分布がマッピングされて良い。ある範囲の古典的記載岩石学的解析が行われて良い。他方図4Bは、3D断層撮像から得られた対応するレジストレーションがなされた薄片を図示している。3D画像における減衰情報を、光学顕微鏡画像から得られた鉱物学上の情報に直接相関付けることに基づいて、顕微手法から得られた鉱物学上の情報は、3DマイクロCT画像データ内部の情報を追加/伝播するように相関付けられて良い。
【００４７】
別な例が図5に図示されている。図5は、50nm解像度で得られたミクロな試料の画像を図示している。図6は、2.5μmの解像度で得られた断層画像からのレジストレーションがなされた薄片を図示している。マイクロCTデータからのグレースケール減衰情報は、50nmスケール情報（たとえば高解像度での孔のサイズ、孔の形状、グレインのサイズ及び形状）に対して直接相関付けられて良い。
【００４８】
3. 当該方法はまた、様々な部位を識別する高解像度の顕微手法を用いる能力をも供する。前記様々な部位は、対応する3DマイクロCT画像と空間的に関連づけられる。一旦これらの部位が、空間的な関係によって、2D画像と3D画像との間での重なり領域について特定されると、残りの3D空間についてもこれらの部位を推定することができる。この手法は有効に顕微鏡データを利用することで、追加の情報を3D画像データに組み込む。図7A-7Cの画像は、どのようにして50nmスケールの2DのSEM画像から得られたグレインサイズと孔のサイズの情報が、図7Aの3D断層画像に埋め込まれた（マッピングされた）孔の構造の情報を生成することができるのかを示している。よって、図7AのマイクロCT画像からのデータと図7B及び図7Cの薄い部分のデータとを結合した合成画像データの組が生成される。そのような合成データは、前記薄い部分内部での微少孔領域の探索によって、ナノメートルサイズの孔のサイズ情報を抽出することを可能にする。nmスケール（50nm以下）で画像情報を探索することが可能であるため、これらのスケールでの3D画像における接続性の特性を測定することが可能となる。
【００４９】
一の応用では、孔のサイズとグレースケール減衰とが直接的に相関することで、本来は単純な有効率として割り当てられていた値を、ボクセルの孔のサイズ情報として3Dデータで伝播させることが可能となる。鉱物学上の情報、音響/弾性特性、及び表面化学情報もまた、最適となるように3DマイクロCTデータを追加する確率論的手法によってマッピングされて良い。
【００５０】
湾曲歪み補正：幾何学的画像補正
顕微鏡の中には、他の画像とのレジストレーションが行われる画像に幾何学収差を導入するものがある。このことはたとえば、画像中における2つの対象物間の測定間隔が、前記2つの対象物間の実際の間隔と一致しないことを意味すると考えられる。正確な定量的情報を推定する前、又は、これらの情報を適切に合成するために2つの画像が正確にレジストレーションされる前に、この幾何学な歪みは除去される必要がある。そのような歪みを画像から除去する処理は湾曲歪み補正と呼ばれる。
【００５１】
湾曲歪み補正は2つの段階からなる。最初に、観察方向に係るワーピング変換が決定される。続いて、逆ワーピング変換が画像に適用されることで、幾何学的に補正された画像が生成される。
【００５２】
ワーピング変換の一例は、正しい相対位置が既知である点からなる組で構成される画像を取得することである。よってこの情報は、歪んでいない画像中の点と歪んでいる画像中の点との間での滑らかなマップを構築するのに用いられる。任意の点の集合がこの目的に用いられて良い。しかし繰り返しになるが、補正の正確さは、その補正の複雑さと滑らかさに依存する。
【００５３】
以降では、2D画像のワーピング変換を自動で決定し得る典型的な方法について記載する。同じ原理は3Dの場合にも適用可能である。
【００５４】
湾曲歪み補正の第1段階は以下を含む。
1. 規則的な正方形グリッドの画像が取得される（図8）。この画像は、グリッド像と呼ばれ、画像中の歪みの決定に用いられる。
2. 前記グリッド像の中心に位置する小さな領域910−各方向において2,3のグリッドセルを含むのに十分な大きさで、かつ実線の暗い影付き四角形で表される−が抽出される。明るい影付き四角形912は、ステンシルパターンが繰り返される場所を示している。
3. 前記グリッド像と前記ステンシルとの間の相互相関が計算され、かつ相関ピークが決定される。前記相関ピークは、同一パターンが再度現れる場所を示している。小さな歪みについては、単純な並進で十分正確である。大きな歪みについては、前記ステンシルの回転及び伸縮が必要になる場合がある。いずれの場合でも、フーリエ変換が、これらの研鑽を効率的に行うのに用いられる。
4. 図10に図示されているように、前記像の中心に最も近い5つの相関ピーク並びに該ピークとグリッドの主軸u及びvとの関係が決定される。このことは、それ自体が中心位置であるだけではなく、前記グリッドの各主軸に沿った第1相関ピークでもある。これらの軸は前記画像の軸と位置合わせされている必要はない。
5. これら5つの相関ピークから、前記グリッドの主軸が決定されて良い。よってこれらの軸は、残りの相関ピークの正しい位置−つまり前記残りの相関ピークが歪んでいない画像中で存在したであろう位置−を推定するのに用いられる。
6. 正しい位置を各相関ピークに割り当てることで、この組についてのワーピング変換が定義される。よって、これらの点の間での位置についてのワーピング変換は適切な補間法によって定義される。
【００５５】
第2段階は前記画像の湾曲歪み補正を行う。この地点では、歪んでいない画像を得るため、前記歪んだ画像に逆ワーピング変換が適用される。つまり、前記歪んでいない画像中の各ピクセルは対応するワーピングされた位置へ移され、かつ前記ピクセルの値は前記歪んでいない画像中での補間によって求められる。前記歪んでいない画像に対して位置合わせされた規則的なグリッドによって処理することで、計算は大幅に単純化される。計算効率のため、この処理は2つのステップとして実行される。
1. 前記歪んでいない画像に対して位置合わせされた規則的グリッド上でワーピングパラメータ（並進、回転等）を計算する。
2. 前記規則的グリッド上での補間によって各像点のワーピングパラメータを求める。湾曲歪みを補正する処理が合成画像に対して及ぼしうる効果が図11に図示されている。図11の観察から、湾曲歪みが補正されたSEM画像から得られた合成画像（図11B）の品質は、湾曲歪みが補正されない場合に得られる合成画像（図11A）の品質よりも優れていることは明らかである。
【００５６】
上述した湾曲歪みの補正方法は例示であり、他の方法を用いることも可能であることに留意して欲しい。
【００５７】
C. 互いにレジストレーションされた複数の3D画像の処理
任意に配向（たとえば回転、伸縮、シフト）する3D画像を3D画像に対して正確にレジストレーションすることに依拠する他の応用についてここで説明する。結果として得られた合成情報は、多くの複雑な実験が有孔性材料について行われた後、ある範囲の条件下で同一の有孔性材料についての物理的特性の測定を可能にすることによって、従来のマイクロCTイメージングの能力を大幅に拡張する。
【００５８】
かなりの長期間にわたって得られた一連の3D画像が3Dでレジストレーション可能である場合、流体相分布は、顕著な劣化の後又は長い平衡時間の後に、3Dにおいて孔のスケールで直接的にマッピングされて良い。よって全範囲の所望の物理特性が画像データに基づいてシミュレーションされて良い。これらには、相対透水率（有孔性材料内部の複数の流体相の透水率）、複数の流体相が集まった岩石の抵抗力、弾性応答特性、NMR応答、及び動的毛管圧力（たとえば非特許文献3を参照のこと）が含まれる。流体分布は、孔ネットワークモデルの結果と直接比較することが可能で、かつ孔スケールのモデル化手法の精度を上げるための校正ツールとして用いられて良い。
【００５９】
開示された方法は概して、孔構造、鉱物相構造、及び3次元有孔性媒質の孔空間内部での複数の流体相の分布の変化の視覚化並びに定量化に関する。有孔性媒体の変化の直接的な定量化から、様々な飽和状態及び様々な表面化学状態下での有孔性材料のある範囲の物理特性が推定されて良い。これらの特性には、様々な流体相の相対透水率、抵抗力、誘電応答、NMR応答、及び弾性率が含まれる。材料には、油層のコア材料、土壌、粒状材料、及び他の材料が含まれる。
【００６０】
開示された方法は、現実的な条件下での多相流、岩石の機構、及び油層岩石物理学上の特性の理解の必須要件を解決する。たとえば、流体の吸水特性はコアの未精製油の経年劣化後に変化することは既知である（たとえば非特許文献7,8を参照のこと）。他の測定もまた、油層温度及び岩圧での実験を行えることが求められる。提案された方法は、説明したアルゴリズムによって実装されたソフトウエアによる、試料の撮像、実験室中での前記試料のある範囲の複雑な実験の実行とそれに続く前記試料の再撮像、及び、前記実験前後で得られた3D画像のレジストレーションを可能にする。
【００６１】
レジストレーションされた複数の3D画像に係るデータの合成方法
この方法の手順が図12に図示されている。当該方法は、関心対象である有孔性材料の3D画像を取得する手順1202で開始される。この撮像は最初の参照状態Aを定義する。一般的には、この状態で研究される試料は、乾燥状態（図13）又は飽和流体（図14）である。画像データは手順1204においてメモリ内に記憶される。
【００６２】
続く手順1206は、孔構造、鉱物相構造、若しくは孔への流体分布を構造的又は科学的に変化させることが可能な有孔性材料上での1つ以上の実験の実行を含む。実験例には以下が含まれて良い。
− 有孔性材料についての流体変位の研究。該研究は以下のようなものであって良い。
・ 非混和性流体の変位
・ 流体の侵襲
・ 界面活性剤、微生物、ポリマー、ゲル、又はコロイドによるフラッディング
・ 混和性流体の変位
− 反応による効果。該効果は以下のようなものであって良い
・ 反応流
・ 溶媒によるフラッディング
・ 堆積及び反応による汚損
・ 生物付着
− 力学的効果。該効果には以下のようなものがある。
・ グレインの鉱物学的特性/破壊
・ 粘土の膨張、付着、及び分離
・ 流体の音響学的シミュレーション
− 孔及び鉱物相の構造への流体の効果
・ 流体の交換
・ 微粒子群、ポリマー、及び/又は水性溶媒によって孔を塞ぐこと
・ 粘土のマイグレーション（天然）
・ コロイドの付着/分離及び汚損
・ 粘土の膨張
− 濡れ性分布。該濡れ性分布には以下が含まれる。
・ 流体分布を孔のサイズ及び構造に相関させること
・ 濡れ性及び孔への流体の分布に対する未精製状態での経年劣化の影響
・ 液体/液体及び液体/固体界面の表面化学特性の測定（たとえば接触角、流体膜の存在）
当該方法は、手順1208において、1つ以上の実験を行った後に前記試料材料を再撮像する手順、及び、手順1204において、前記の再撮像された画像を記憶する手順をさらに有する。前記試料は複数回再撮像されて良い。前記試料による実験を行う代わりに、手順1206は、より高い解像度での試料のさらに小さく分割された部分の取得に係るものでも良い。その結果、実験前後で取得された画像又は様々な撮像解像度で得られた画像のいずれかを表す1対の画像が生成される。よって2組のデータは、3D画像のボクセルレジストレーションを実行することによって得られた画像のうちの2つから得られた3Dデータを位置合わせして、かつ重ねることによって、合成画像情報を得るのに用いられて良い。そのような位置合わせ及び合成は、続いて別な画像の対の間で行われても良い。2つ以上の画像の情報を合成する処理は、両画像からの情報に基づいて、前記試料の特定領域の物理特性又は化学特性を決定する手順を有して良い。当該処理はまた、孔の構造、鉱物相構造、又は孔への流体分布の変化を定量化する手順をも有して良い。
【００６３】
ここで図13と図14を参照すると、図は、流体取り込み後の有孔性岩における空気/水相を図示している。3Dレジストレーションが実行されることで、流体の孔スケール分布を3Dで直接的に視覚化及び定量化することが可能となる。特に図14は、乾燥した岩へ濡れ性流体が取り込まれる実験を行った後、図13から得られた岩石の流体分布を図示している。
【００６４】
図15は、一対の3D画像データファイルの画像レジストレーションから得られた流体分布を表す3Dビデオからのスナップショットを図示している。
【００６５】
試料材料の再撮像は、各実験が行われた後又は一連の実験が行われた後、何度も繰り返されて良い。繰り返しになるが画像データは、コンピュータのメモリ又はデータ処理を取り扱う他の電子システム（通常はハードドライブ）に記憶される。よって、異なる画像を取得する間の期間に制限はない。よって2つ以上の記憶された画像のボクセル−ボクセルレジストレーションは、任意の数の画像から得られる情報を結合させることで、全実験値に追加することを可能にする。
【００６６】
任意の期間にわたって行われる任意の数の実験から合成画像を生成するために提案された方法は、従来技術を有する2つの実験方法とは異なる。非特許文献9に代表される前者の方法では、全ての実験はビームライン中で行われる。このことはつまり、実験の性質及び期間は、ビームラインの利用可能性、空間についての考慮、並びに、温度及び機械的安定性の要求によって制限される。後者の方法（たとえば非特許文献10を参照のこと）では、試料が撮像され、続いて実験の修正を行うために移され、その後再撮像のためにビームラインに戻される。しかし異なる画像のレジストレーションを行うための努力は行われない。この場合では、合成画像の情報は直接得られず、唯一の結果は空間的な平均である。
【００６７】
最終手順1212は、手順1210においてレジストレーションされた3D画像から得られた合成情報を用いることで、孔の構造、鉱物相構造、又は孔への流体の分布の変化を定量化する手順を有する。これらの変化の有孔性材料の物理特性に対する示唆が調査されても良い。たとえば濡れ性分布の孔スケール分布の直接的な定量化から、様々な飽和状態及び様々な表面化学条件下でのこれらの有孔性材料のある範囲の物理特性の直接的なシミュレーションが検査されて良い。関心特性には典型的には、様々な流体相の相対透水性、抵抗力、誘電応答、及び弾性が含まれる。
【００６８】
図16A及び図16Bは、孔のサイズ分布と空気/水流体の含有量との間の相関を示している。3Dレジストレーションにより、孔スケールでの流体の分布を3Dにおいて直接的に可視化及び定量化することが可能となる（図16A）。図16Bは、濡れ性を有しない（空気）相が3Dにおいて存在する場所の分布を示している。この場合では、系は本来水で濡れているので、残りの空気相は基本的に大きな穴にとどまり、かつ水相は小さな孔に集中する。
【００６９】
D. 典型的な応用
孔の構造、鉱物相の構造、若しくは孔への流体の分布の起こりうる構造変化又は化学変化を定量化することが可能であることが重要となる、ある範囲の学問分野にわたる有孔性材料に関連する多数の実験がここで論じられる。
【００７０】
1. 孔のスケールでの有孔性材料の流体変位の記述。有孔性材料の撮像の第1用途の中には、様々な濡れ性条件及び飽和状態の下で記載され、かつ再撮像されるものがある。複数の相が有孔性材料内部で撮像されて良い。流体変位特性に関連する実験には以下が含まれて良い。
【００７１】
a. コア解析における非混和性流体変位（上流の油とガス/地下水及び水門学）
・ 複数の飽和状態での排水及び吸水フラッディング
・ 様々な濡れ性状態（油が濡れる、水が濡れる、混合物の濡れが大きい、混合物の濡れが大小さい）下での排水及び吸水
・ たとえばフラッディング中及びフラッディング後における、未精製油の経年劣化が、炭化水素、ブライン、ガスなどに及ぼす効果の探索
・ ある範囲の流体分布の探索。
− 油/炭化水素
− 水/ブライン
− CO2
− 泥
− 入り込む流体
− ポリマー
− コロイド
− 分散剤
− 孔へ入り込む流体（たとえばポリマー、ゲル、微生物種）
・ 意図した充填の試験
− 孔構造における抽出を最大にするための最適なブライン
− 様々な油相/ブラインの混合物での孔スケール分布の試験
− 改良/改善された油の回復戦略（たとえばWAG注入及びコロイドシミュレーション）の効率
・ コア材料への泥の侵襲（損傷の生成）の撮像
b. 混和性流体変位
・ 混和性フラッディング(EOR)
・ 地下水の浄化
− 孔スケールでの毒物又は放射性物質の広がりの分布の定量化
− 汚染物（化学的、水門学的）の試験方法
2. 孔及び鉱物相の構造への流体の効果：孔の構造、鉱物の構造、及び有孔性材料の濡れ性に対する流体の交換の効果（たとえば孔の幾何学形状に対するブラインの希釈効果）
a) 粘土のマイグレーション（天然）
b) コロイドの付着/分離及び汚損
c) 粘土の膨張
3. 濡れ性の分布。多数の流体相によって占められる有孔性材料の流体/流体界面及び流体/固体界面の濡れ性状態が評価される。これらには以下が含まれる。
a. 流体分布を孔のサイズ及び構造に相関させる
b. 液体/液体界面及び液体/固体界面の表面化学特性の測定（たとえば接触角、流体膜の存在）
4. 反応の効果。有孔性材料の孔及び鉱物相の構造への反応種の効果を直接孔のスケールで調査する。流体の濡れ性についての研究は複数回の手順で行われても良い。
a. 反応流
− 超臨界CO2
− 酸のフラッディング
− 溶媒に基づくフラッディング
b. 堆積及び反応による汚損
c. 生物汚損
5. 力学的効果。有孔性材料の孔及び鉱物相の構造への力学的な応力と歪みの効果。これらには以下が含まれる。
− グレイン鉱物学的効果/破壊
− 粘土の膨張、付着、及び分離
− 流体の音響試験。炭化水素の回復を改善するための音響波の試験的利用。
【００７２】
E. 画像レジストレーションのワークフロー
以降で説明するレジストレーションのワークフローは、2D-3D画像レジストレーションにも3D-3D画像レジストレーションにも適用可能であることに留意して欲しい。定義を多少変更させれば同一のレジストレーションのワークフローは2D-2D画像レジストレーションにも適用可能である。
【００７３】
定義
Tを1つの独立の画像全体を表すものとし、Ti,j,kを指数(i,j,k)での前記独立の画像の強度を表すものとする。Fを他の独立の画像全体を表すものとし、同様にFi,j,kを指数(i,j,k)での前記他の独立の画像の強度を表すものとする。さらに、I(n)を、dn倍のダウンサンプリングを行うことによって独立の画像Iから生成される独立の画像であるとする。T(0)⇔T及びF(0)⇔Fとなるように常にd0=0と定義する。独立の画像Tは動く独立の画像で、かつFは固定された独立の画像である。Tが動くことで、Tは固定された独立の画像に対してレジストレーションされる。たとえばTは、3DマイクロCT画像であるFに対してレジストレーションされる2DのSEM画像であって良い。Tが2D画像である場合、Tは3次元空間内の一面であると単純に解される。つまり第3の指数は、k=0に固定される。あるいは別な例として、Tは、3Dの乾燥状態のマイクロCT画像であるFに対してレジストレーションされる3Dの湿った状態のマイクロCT画像であって良い。独立の画像FとTはまた一対の2D画像であっても良い（たとえば図17に図示されたSEMの薄い部分の画像の一部分174）。またわずかな回転によって変化させることで、後述の方法は2D-2Dレジストレーションにも適用される。
【００７４】
Ωt⊂R3であるとき、t(n):Ωt→Rを独立の画像T(n)の補間された画像T(n)であるとする。同様に、Ωf⊂R3であるとき、f(n):Ωf→Rを独立の画像T(n)の補間された画像F(n)であるとする。
【００７５】
Ψ(n):Rm(n)×R3→R3を空間変換と定義する。Ψ(n)は、一組の変換パラメータφ∈Rm(n)が与えられると、3Dの点を3Dの点へ写像変換する。
【００７６】
M(n):Rm(n)→Rを、2つの補間された画像t(n)とf(n)との間でのレジストレーションの品質を示す距離と定義する。またM(n)の値が小さければ小さいほどレジストレーションは良好で、M(n)の値が大きければ大きいほどレジストレーションは不十分である。
【００７７】
初期の独立の画像t(0)とf(0)が与えられると、2つの画像をレジストレーションさせる変換パラメータφを求める必要がある。これは最適化問題に再定式化することができる。
【００７８】
レジストレーション最適化問題：空間変換Ψ(0)が与えられると、距離M(0)及び画像t(0)とf(0)によって、φ*∈Rm(0)が求められる。このとき、
【数１】
が満足される。
【００７９】
一般に、M(0)は多くの極小値を有する複雑な関数で、かつ（特にフル解像度の画像について）評価するのには計算上手間がかかる。従って、局所的な数値最適化法（シンプレックス法、パウエル法、勾配降下法等）単独ではφ*（の近似値）を巧く求めることができない。以降では、レジストレーション最適化問題を解決することが可能な多解像度多出発点の大域最適化法が与えられる。
【００８０】
一般的方法
1. ダウンサンプリング因子d0,…dN、空間変換Ψ(0),…Ψ(N)、距離M(0),…M(N)、及び最善の変換パラメータの数J(0),…J(N)を選ぶ。たとえば以降のこれらの変数の選択を参照のこと。
2. φ(N)⊂Rm(N)を、Rm(N)から得られた有限個のサンプリング点と定義し、かつP(N)を、組φ(N)中の要素数とする。各φP(N)∈Φ(N)、p=1,…P(N)について、mp(N)=M(N)(φp(N))を評価する。以降のたとえばΦ(N)を参照のこと。
3. 組Φ(N)（ハット）={φj:φj=φpj(N)、j=1,…J(N)}を生成する。ここでpjは、mpj(N)がすべてのmp(N)、p=1,…P(N)のうちでj番目に小さな値であるように定義される。
4. n=Nと設定する。
5. n=n-1と設定する。n≧0の場合には6へ続き、そうでない場合には9へ進む。
6. j=1,…J(n+1)の各々について、開始点φj(n+1)∈Φ(n+1)（ハット）を有する局所的数値最適化アルゴリズム（たとえばパウエル法）を用いてM(n)の極小値を計算する。極小値を計算するための変換パラメータをφj(n)（ハット）とし、かつこれらの極小値での距離をmj(n)（ハット）=M(n)(φj(n)（ハット）)とする。
7. 変換パラメータの組Φ(n)（ハット）={φj:φj=φpj(n)（ハット）、j=1,…J(n)}を生成する。ここで、pjは、mpj(N)（ハット）がすべてのmp(N)（ハット）、p=1,…P(N)のうちでj番目に小さな値であるように定義される。
8. 5へ進む。
9. 次式を、t(0)の画像をf(0)の画像へレジストレーションするのに最適なパラメータの組と設定する。
【００８１】
【数２】
10. 前記最適なφ（ハット）の変換パラメータを用いることによって、レジストレーションされた状態である独立の画像の新たな対を生成する。
【００８２】
手順1は、ダウンサンプリング因子、許容される空間変換、及び所与の組の変換パラメータで実現されるレジストレーションの品質評価に用いられる距離の種類についての詳細を定義する。典型的には、ダウンサンプリング因子は2の乗数となるように選ばれる。なおd0=1で、n0<n1であればdn0≦dn1である。例の結果、d0=1,d1=2, d2=4,d3=8,
d4=16,及びd5=16である。最大のダウンサンプリング因子は、特徴部が最低解像度の画像となるのが防止されるように選ばれなければならない。
【００８３】
空間変換Ψ(0),Ψ(1),…Ψ(N)は、どのようにしてt(n)の領域Ωtからの座標が、f(n)の領域Ωfへレジストレーション/重ね合わせされるのかを定義する。一の可能性は、すべての空間変換が等しくなるように、Ψ(0)⇔Ψ(1)⇔・・・⇔Ψ(N)⇔Ψ（オーバーバー）と選ぶことである。これは、図中の例をレジストレーションする場合には正しく、かつΨ（オーバーバー）は3D相似変換となるように選ばれた（Ψ（オーバーバー）:R7→R、7の自由度は、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータを有する）。
【００８４】
選択可能な実際の距離の指標には様々なものが存在し、繰り返しになるが、妥当な選択とは、各探索レベルn,M(n)（オーバーバー）の距離について同一の数式を有することである。つまり探索レベル間での唯一の違いは、距離が対応するダウンサンプリングされた画像対を用いて評価されることである。有用なコストの指標には、相関係数、相関比、及び（規格化された）相互の情報が含まれる。図中のレジストレーションされた例については、M(n)（オーバーバー）が負の相関係数として選ばれた、
【数３】
ここでt(n)（オーバーバー）は平均で、σt(n)は試料座標xlでのt(n)画像の強度の標準偏差で、f(n)は平均で、かつσf(n)は試料座標Ψ(φ,xl)でのf(n)画像の強度の標準偏差であり、試料の空間座標は画像領域の交差部分xl∈(Ωt∩(Ψ-1(φ)（オーバーバー）・Ωf))、から取得される。ここでΨ-1（オーバーバー）はΨ（オーバーバー）の逆関数で、Ψ（オーバーバー）（Ψ-1（オーバーバー）(x)）=x,∀x∈R3となる。
【００８５】
J(n)の値は、探索レベルnから探索レベルn-1へ渡される最善のランク付けがなされた変換パラメータの数である。典型的には、J(0)=1で、かつn0<n1のときにはJ(n0)≦J(n1)となる。J(n)の結果の一例は、J(0)=1, J(1)=1, J(2)=1,
J(3)=4, J(4)=32, J(5)=64である。
【００８６】
手順2では、目標は、所望の変換パラメータφ*の組に「近い」変換パラメータの組をすくなくとも1組求めることである。この手順では、包括的な型の探索が、より高い解像度での所望のレジストレーションを実現すると思われる変換パラメータの組を発見するために実行されて良い。典型的には、Φ(N)の組は、変換パラメータΦ(N)={φα:Zm(n),α=0,…A}の正方形グリッドとして定義される。ここでα=[α0,α1,…αm(N)]Tは多重添字である。φαは次式に示す行列で表される。
【００８７】
【数４】
ここでφmin∈Rm(N)は変換パラメータ成分の下限で、δ∈Rm(N)は変換パラメータ空間内のグリッドのステップサイズである。最小の並進範囲は、Ωf∩Ψ(N)(φ)・Ωtが全てのφ∈Φ(N)に対して空集合とならないことを保証するように選ばれて良い。ここで典型的な並進ステップサイズはダウンサンプリングされたボクセルサイズである。回転値は-180°から180°の範囲で、ステップサイズは約3°である。手順2は基本的には、距離の大域最小の「捕獲半径」の範囲内にある1組の変換パラメータを求めるために低解像度の画像に適用される大域最適化法である。この段落は、「包括的な探索」による大域最適化法について説明している。しかし最適レジストレーションパラメータに「近い」変換パラメータ（又はその複数の組）を発見するため、他の大域法も手順2において適用されて良い。
【００８８】
手順3は、最高のレジストレーションから最低のレジストレーションまでの変換パラメータの組を順序づける。繰り返しにより、手順4から手順8は、順次より高解像度の画像へ局所最適化を実行することによって、手順2で求められた変換パラメータの組を徐々に改善する。
【００８９】
手順9は最適なレジストレーションを与える変換パラメータを選ぶ。手順10は、画像の再サンプリングを（その画像の最適変換パラメータを用いて）行うことで、後続の解析を行うために直接比較することが可能な一対のレジストレーションされた独立の画像を生成する。
【００９０】
アルゴリズム
上述のレジストレーション法を実装するためのアルゴリズムが図17に図示されている。図17では、3D画像1702と2D画像1708との間でのレジストレーションが実現される。最初に得られた画像1702は、フィルタリング、マスキング等によって、プロセッサ3605Aによって前処理されることで、レジストレーションの準備が整った画像1704が得られる。同様に、個々の画像1706もまた前処理されることで、2Dの操作された画像が得られる。特に、2Dの再分割された部分1706の前処理は、湾曲歪みの補正及び分割画像を一つにまとめることで1つの画像1710を得る手順を有して良い。
【００９１】
図17のアルゴリズムの手順1では、2D画像1710を3D画像1704に対して位置合わせするため、相似変換がプロセッサ3605Aと3605Bによって選ばれる。7の自由度にはこの変換も含まれる。具体的には、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータである。負の値を有する相関係数は、位置合わせの品質を評価するのに用いられる。-1の値は完全な相関を意味する。
【００９２】
手順2では、一対の低解像度の画像が、プロセッサ3605Aと3605Bによって用いられることで、画像1704と画像1705との間での位置合わせについての一組の初期推定が得られる。手順2の実装は、タスクの並列化によってその手順の実行時間を減少させるという利点を得ることができる。この用途におけるタスクの並列化の実装についてのさらなる詳細は本明細書においてさらに与えられる。手順3では、局所最適化用の開始点として用いられる一群の最善の変換パラメータが選ばれる。当該処理のこの段階では、画像1704と画像1710とは最低解像度の状態で大まかな位置合わせがなされている。
【００９３】
手順4では、多数の探索に関してのループカウンタが、プロセッサ3605Aと3605Bによって初期化される。続く手順5では、ループカウンタは、プロセッサ3605Aと3605Bによってすぐに1つずつ減ってゆく。ループカウンタが0以上の値を保持する場合、当該方法は手順6へ進む。手順6では、プロセッサ3605Aと3605Bは、現時点での高解像度画像で局所最適化を用いることによって、過去の探索レベルn+1からの位置合わせを精緻化する。手順6の実装は、低解像度の画像についてのタスクの並列化を行うという利点を得ることができる。前記タスクの並列化では、3605Aと3605Bがそれぞれ異なる複数の開始点の部分集合を処理する。より高解像度の画像の対については、画像に必要なRAMが1つのコンピュータプロセッサの物理RAMを超える場合、データ並列化を実装することができる。前記データ並列化では、3605Aと3605Bはそれぞれ、各プロセッサからアクセス可能なメモリ内に記憶された各プロセッサに対応する部分像に関する指標の各独立した一部分を計算する。係る実装についてのさらなる詳細は本明細書においてさらに供される。手順7では、プロセッサ3605Aは、次の探索レベルの開始点として用いられる現在の探索レベルから最善の変換パラメータの群を選ぶ。手順8では、そのループは手順5へ戻ることによって続けられる。
【００９４】
nが負の場合には他のルートがとられる。この場合では、手順5の後に、プロセッサ3605Aと3605Bは手順9へ進む。そこでプロセッサ3605Aと3605Bは探索レベルn=0から最善の変換パラメータを選ぶ。これらは、フル解像度での最低指標値を供するパラメータである。これらのパラメータを用いることで、手順10では、プロセッサ3605Aと3605Bは、3D断層画像を再サンプリングすることで、2D画像に対応する薄片を生成する。この結果、同一サイズである2つの位置合わせされた2D画像が生成される。よってこれらの画像は、後続の定性的解析及び/又は定量的解析を行うのに直接用いられて良い。
【００９５】
図17の手順1から手順10は、図1で開示された画像処理方法の高次表現に係る手順108と指称されるレジストレーション方法を実質的に規定する。同様に図18は、図12で開示された画像処理方法の高次のダイアグラムに係る手順1210を示している。
【００９６】
図18に示されたレジストレーション方法は、2つの3D画像1802と1808にも適用可能である。画像1802は乾燥状態のコア試料の3DマイクロCT画像であって良い。対照的に、画像1808は、特別の実験を行った後の同一試料の3DマイクロCT画像であって良い。前記特別の実験にはたとえば、特別の流体−たとえば水、石油、又はガス−へコアを導入する手順が含まれて良い。画像1802と1808のいずれも、たとえば単純な閾値設定を利用したフィルタリング及びマスキングによって前処理される。前記前処理の結果、画像1804と1810がそれぞれ得られる。繰り返しになるが、手順1では、プロセッサ3605Aは、2つの断層画像1804と1810の位置合わせを行うため、相似変換を用いる。繰り返しになるが7の自由度がプロセッサ3605Aによって選ばれる。前記7の自由度は、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータを有する。負の重なり割合が用いられる。この負の重なり割合は、両画像での同一のマスキングされていない強度を有する座標の割合を単純に計数するものである。-100%の値は両画像においてマスキングされていない領域が完全に重なっていることを意味する。
【００９７】
手順2から手順9は概して図17に図示された手順と同一である。繰り返しになるが手順10は、最終手順であり、かつ、乾燥状態の3D断層画像1804に対応する領域を生成するため、湿った状態の3D断層画像1810を再サンプリングするプロセッサ3605Aと3605Bを有する。その結果、一対の3D画像1814と1816が、プロセッサ3605Aと3605Bによって得られる。一対の3D画像1814と1816は、同一サイズであって、かつ定量的な解析目的のために合成されてよい。
【００９８】
図19は、図17及び図18の手順2の様々な下位手順をより詳細に図示している。手順2は、7次元の変換パラメータ空間において規則的なグリッドをプロセッサ3605Aと3605B によって生成する手順を有する手順2aで開始される。典型的には並進パラメータのステップサイズは最低解像度の画像のボクセルサイズ長である。z軸の回転パラメータのステップサイズは2°である。x軸とy軸の回転パラメータ及び伸縮パラメータは（それぞれ0°、0°、及び1で）固定されて良い。一旦手順2aから規則的なグリッドが生成されると、プロセッサ3605Aは、kを7Dパラメータグリッドにおける点の数とすることによって、パラメータの初期化を行う手順2bへ進む。続く手順2cでは、パラメータkは1だけ減少し、かつプロセッサ3605Aと3605Bによって、新たなパラメータ値が0と比較される。
【００９９】
減少したkが0以上である場合、プロセッサ3605Aと3605Bは、手順2cから手順2dへ進む。手順2dでは、プロセッサ3605Aと3605Bはそれぞれ、7Dグリッド点の数kの変換パラメータの指標を評価する。前記指標の値がこれまでに得られた最善の指標値よりも低い場合、プロセッサ3605Aと3605Bは、最低の指標値を与える変換パラメータを変換パラメータ数kに置換する。手順2dの結果は、プロセッサ3605Aと3605Bを手順2cへ戻すためである。あるいはその代わりに、パラメータkの値が、手順2cで減少したように、負である場合、プロセッサ3605Aと3605Bは手順2cから手順2eへ進む。手順2eは手順2の全体を終了させる。パラメータkが負の値をとるということは、7Dグリッドの包括的な探索が完了し、かつ指標となる大域最小値の捕獲半径の範囲内に属する少なくとも一組の変換パラメータが特定されたことを意味する。ここでタスクの並列化は、各異なるプロセッサ−たとえば3605A、3605B等−間のkのグリッド点を分割することによって実装されて良い。各プロセッサは、グリッド点の部分集合についての指標を評価する。最善のパラメータの群は、各プロセッサから収集され、かつ照合されることで、最善の変換パラメータからなる最終群が生成される。
【０１００】
並列実装
上述の大域最小化法を計算機上で実行可能にするため、高性能非均質メモリアクセス(NUMA)プラットフォーム上のソフトウエアで上述の手順を実行されるときに使用される2つの並列化法が存在する。第1の方法はタスクの並列化法である。前記タスクの並列化法では、個々のタスクは計算ユニット上で独立に実行される。第2の方法はデータの並列化法である。前記データの並列化法では、各計算ユニットは独立の画像対データの部分集合だけしか含んでいない。
【０１０１】
タスクの並列化法は、上述のワークフローの手順2において距離の値mk(N)を計算するときに大いなる利点となる。低温では、各計算ユニットは、他の計算ユニットとは独立してmk(N)の値の部分集合を評価することができる。全ての計算ユニットが、それらのmk(N)の値の部分集合を評価したとき、結果として収集段階では、手順3で実行されるように、最善の変換パラメータを順位付けするため、各計算ユニットの各々からmk(N)の値の部分集合が集められる。同様に、各計算ユニットが極小値Φj(n)（チルダ）を計算するため、他の計算ユニットとは独立にして繰り返し最適化の一部を実行することができる手順6で利点を得るためにタスクの並列化が用いられて良い。繰り返しになるが、結果として実現される収集段階は、手順7で実行されるように、各計算ユニットからのmk(N)の値の部分集合を収集して順位付けするのに用いられる。タスク並列化は低解像度の画像（たとえばt(N), t(N-1), f(N),及びf(N-1)）にとって有利である。なぜなら基礎となる独立の画像対の全データは計算ユニットの局所RAM内に記憶することができるからである。
【０１０２】
高解像度の画像（たとえばt(0), t(1), f(0),及びf(1)）にとっては、データの並列化法が好ましい。この手法では、各計算ユニットは独立の画像のデータの部分集合しか有していない。よって手順7は、マスターワーカー並列化によって実行される。マスター計算ユニットは、繰り返し最適化アルゴリズムを制御する。前記繰り返し最適化は、極小値を決定するため、指標となる距離の多重評価を必要とする。指標となる距離の評価は、分配して、かつ並列に行われる。各ワーカー計算ユニットはまた、該ワーカー計算ユニットの局所RAM内に保持された画像データの部分集合に基づいて、全体の指標となる距離の一部を計算することができる。
【０１０３】
E. さらなる応用例
X線マイクロ断層撮像(CT)は、ミクロンスケールでの孔の構造の詳細な3D画像を生成する能力を有する。X線断層撮像は、検出されたX線の減衰を介して位相情報を供する。このデータについては詳述しないが、様々な2D顕微鏡手法によって供されるものとして有用である。たとえば光学顕微鏡及び走査型電子顕微鏡は、試料の詳細な鉱物学的内容の情報を供することができる。走査型音響顕微鏡は、鉱物の音響（弾性）特性のマッピングを可能にする。2D表面又は薄い部分での実行が可能な他の顕微鏡手法は、研究対象であるコア材料（試料）の材料特性の理解にとって重要な他の情報（たとえば鉱物学的知見、表面特性、局所的表面相互作用等）を与えることができる。一の欠点は、得られたデータが2Dであることである。2D方法によって利用可能な情報が3Dにおいて改善された情報として得られることは、大きな利点である。ここで述べている方法は、2D画像から得られた情報を3D断層画像データの組に定量的に伝播することを可能にする。
【０１０４】
上述のレジストレーション方法はこの伝播処理の一部である。図20は2D画像2010を3D画像2012へレジストレーションした結果を図示している。特に、当該方法において定義されている全ての並進、回転、ワープ、及び伸縮変換が行われた後、2D顕微鏡で得られた薄い部分2010は試料2012内部でレジストレーションされることで、CT画像2012から得られた3Dデータ組への、2D画像2010から得られたSEMデータの結合が可能となる。
【０１０５】
マイクロ有孔性岩石の試料の3D断層画像への高解像度2DのSEM画像から得られた画像データのレジストレーションの例が図21に図示されている。2DのSEM画像2112から得られたマイクロ有孔性データを3D断層画像2110へ結合することによって得られた合成データは、2DのSEM試料での解像度、及び3D断層撮像によって可能な解像度よりも高い解像度で、材料の部位及び特性を識別することを可能にする。一般的には、断層撮像手法は、本質的な解像度の限界（たとえば1μm）を有する。顕微鏡はnmスケールでの探索を可能にする。よって、画像レジストレーションによって2つの画像を合成することで、最大1000倍の解像度のオーダーが利用可能となる。それに加えて、2つの手法の解像度の大きな差異は、3D画像の真の解像度を測定するのに用いることができる。
【０１０６】
図22は、岩石試料の3D断層画像から得られた2D薄片2210と、SEMから得られた2D薄片との直接比較を図示している。断層画像2210は2.1μmのボクセルサイズのX線CT撮像によって得られているのに対し、SEM薄片2212は0.25μmの解像度で得られている。2つの画像2210と2212との定量的な一致が図23に図示されている。具体的には、図23は、高解像度の2D画像2212から得られた（特定の領域の一部に係る）有孔率データを3D薄片2210から得られた減衰データへマッピングするグラフを図示している。グラフ2314から分かるように、減衰データと有効率データとの間には良好な相関が得られた。このため、高解像度2D顕微鏡像から得られた情報（この場合は有孔率）を、3D全体へ伝播させることが可能となる。これは、3D画像の減衰データに基づくものであり、かつグラフ2314によって図示された直接的かつ定量的にマッピングされた相関を用いている。
【０１０７】
2DのSEM画像から高解像度で得られる他の特性もまた、レジストレーションされた3D画像と定量的に合成されて良い。図24A及び24Bに図示された例では、SEM画像2410は0.25μmのピクセル解像度で得られる。画像2414は、2.5μmのボクセルサイズで得られたX線CT画像からのレジストレーションされた薄片である。この画像でのX線減衰は、画像2410内の孔のサイズデータに相関付けられて良い。画像2412は、画像2410から得られた画像データ上での局所的な孔のサイズ分布のマッピングを表す。画像2412中でのグレースケールは、2410から得られた孔のサイズに対応する。明るい領域は大きな項のサイズに対応する。
【０１０８】
図24Bは、局所的な孔のサイズ分布（y軸）とX線減衰からのマイクロ有孔率マッピング（x軸）との相関のプロットを図示している。この相関により、3D画像全体への孔のサイズをマッピングすることが可能となる。
【０１０９】
図25は、2D顕微鏡画像から3D断層画像への定量的な画像レジストレーションの別例を図示している。この例では、元素マッピング（画像2512）から得られた鉱物相の情報が、画像2510によって表される3D画像から得られたレジストレーションされた薄片で得られた減衰データに関連づけられている。2D画像2512に示された詳細な鉱物相の分布は、岩石試料内でのドロマイト及び硬石膏の存在を示している。これらの相は、3D画像2510では明瞭に差異化するのが難しい。両画像間でのレジストレーションを利用することによって、定量的マッピングが、3D画像全体での鉱物相の情報を識別する能力を改善させた。
【０１１０】
ここで、3D-3Dレジストレーションに係るさらなる応用について説明する。係るレジストレーションによって、1つ以上の実験が試料上で行われた後に、その試料中に生じた変化を定量的に追跡することが可能となる。3D-3Dレジストレーションを含む典型的な応用は通常、試料−特別な状態Aにある−の撮像で開始される。この最初の画像が得られた後、試料は撮像装置から取り外され、かつ、様々な実験−たとえば圧密、圧縮、分解、流体変位−が様々な条件下で行われる。係る実験によって試料の状態（たとえばマイクロ構造、鉱物学的状態、流体飽和）が変化する。続いて、この時点で新たな状態Bにある試料は、撮像システムへ再導入されて再撮像される。続いて新たに得られた画像が先に得られた画像に対してレジストレーションされる。2つの画像から得られたデータは合成され、生じた変化は定量化することができる。
【０１１１】
以下の例は、一組の実験中での分解と粒子運動の視覚化及び定量化を有する方法の応用、並びに様々な化学条件下での複数の流体飽和状態の撮像を表している。また複数のスケールでの3D画像の位置合わせ例も与えられている。
【０１１２】
図26及び図27は、一連の実験中に試料体積内部の物理構造において生じた変化を定量的に追跡する実験を表している。係る塩化は鉱物に固有であって良い。岩石試料2610は、実験室で受け取られたような本来の状態である。薄片2616は、試料2610が本来の状態であるときに、その試料2610の3D画像内部から得られたものである。続いて試料2610には、中間状態2612によって示されているように、1つ以上の実験が行われる。図26の場合、試料2610は特定の組成を有する流体に浸される。その流体は、特定の期間中高温かつ高圧に維持される。特別なセル内で維持されたそのような条件は、試料材料の分解及び沈降を起こす。その結果、長期間にわたって試料のミクロ構造が変化する。そのような変化を評価するため、試料は新たな状態2614で再撮像される。3D画像内部からの薄片2618もまた図26に図示されている。2616及び2618のより詳細な拡大図が、2720と2722にそれぞれ図示されている。薄片2720と2722とを比較すると、試料の形態において視覚上明確な差異が存在する。前記差異は特定の実験条件によって引き起こされたものである。薄片2720と2722に係る2つの大きな3D画像を合成する結果、多数のデータが、遷移期間中に試料体積全体にわたって生じたミクロな変化を捉えることができた。例の小さな部分集合は、2616から2618への遷移中に発生した、孔構造、グレイン構造、又は、流体透水率の変化の列挙を有して良い。これらの局所スケールの差異/変化は、当該方法に含まれる直接レジストレーションによってしか識別できない。
【０１１３】
図28は、カルボン酸への曝露後における岩石の分解/沈降特性の変化を研究するための様々な実験を図示している。薄片2810は初期状態における岩石試料の本来の3D画像を表している。薄片2812は、長期間にわたってカルボン酸へ曝露した後の同一岩石試料からとったレジストレーションされた3D画像から得られた薄片を表す。前記曝露の前後で2つの3D画像のレジストレーションを行うことによって、岩石試料中で生じた形態変化を識別することができる。この特別な例における変化は顕著でかつ明らかである。たとえば、本来の画像2810において視認可能な部位2814は、酸への曝露後に撮られた画像から得られた、同一で定量的にレジストレーションされた薄片2812から消滅した。よって、2つの画像間での差異を示す合成された3Dデータの組は、材料構造の変化を定量化することを可能にする（たとえば孔の構造、孔の形状、グレインの形状、ミクロ有孔率）。物理的特性の変化もまた定量化することができる。たとえば分解及び再沈降による試料の流動特性の変化は、合成された3Dデータ組を用いることによって探索可能である。地質材料の分解を調査するための様々な方法について検討するため、さらなる実験が様々な条件下で行われても良い。
【０１１４】
図29は、オイル、ガス、及び/又は水をコア試料全体へフラッディング（循環）させる間に生じる様々な状態によって生成される一連の曲線2910を図示している。各曲線2910は単一の濡れ性状態についての毛管圧力を示す走査曲線である。画像レジストレーションによって、試料のフラッディングサイクル（たとえば1次排出、1次自発吸収、1次強制吸収、2次排出、2次吸収などの後の流体飽和状態）中の流体飽和状態において生じる変化を3Dで視覚化及び定量化することが可能となる。グラフ2910は、特定の一連の実験中に起こる様々な考えられ得る経路−たとえば経路A→B→C→D−を識別している。この一連の実験でわずかな変化が起こった結果、経路がA→E→F→Gとなることもあり得る。一組の実験は以下の手順のうちの少なくとも1つを有して良い。
1. グラフ2910の状態Aでは、岩石試料は水によって十分飽和している。
2. 1に続いて、前記岩石へのオイル又はガスの注入が、特定の最終注入圧力になるまで行われる。その結果、たとえばB（高圧）又はE（低圧）のような状態となる。
3. その後、系を濡れる流体と接触させながらその系を緩和させることによって、前記オイル又はガスは「自発的に」回復する。それによってC又はFのような状態がそれぞれ得られる。
4. その後前記試料は、強制的な濡れ性流体の注入を受ける。その結果、状態D又は状態Dになることができる。
5. 当該サイクルは手順2から任意の回数繰り返されて良い。
【０１１５】
各サイクルの各段階での流体構造の差異を定量化することができる。それにより調査対象である多相流動特性を推定することが可能となる。上述したように、2910は単一の濡れ性状態について毛管圧力を示す走査曲線を示している。複数のサイクルの実験が、他の流体によって、かつ様々な条件下において、同一のコア材料上でさらに行われて良い。そのような実験が行われることで、異なる一組の毛管圧力/水飽和の走査曲線が得られる。これらの状態は撮像することが可能である。流体飽和状態の差異が撮像され、かつ定量的に比較することができる。図30-31に図示された飽和状態は、毛管圧力走査曲線2910上において示された点A-Dに基づく。図30は、状態Aと状態Bとの間での流体分布の変化を図示している。状態Aでは、試料は水3010で飽和している。状態Bは、圧力をかけた状態で、遺留水3012が飽和するまえオイル（又はガス）を試料へ注入した後のレジストレーションされた3D体積から得られた、同一の薄片である。3012では、オイルが試料へ入り込む孔が黒で、かつ残りの水の飽和は小さな（わずかに明るい）孔である。定量的なレジストレーションは、状態Bで注入された流体の量を特定する−これは体積測定に基づいて測定することができる−だけではなく、両流体の詳細な孔スケール分布を直接観察することも可能にする。2つの流体の相対流動能力は、各相における相対透水率のシミュレーションによって3D画像で定量化することができる。そのような定量的な解析は、岩石の初期の炭化水素のフローポテンシャルの推定を行うことができる。前記フローポテンシャルは重要な製造上の示唆となるパラメータであると考えられる。
【０１１６】
図31及び図32Aの状態CとDは、状態Bで観察される初期（遺留水）の飽和からのオイル又はガスの自発変位及び強制変位後にそれぞれ明確になる。ここでは、「自発的な」流体の注入と強制的な水の注入の両方によってオイルが回復した後、どの程度オイル又はガスが試料中に残っているのかを観察することが可能である。この観察は、残りのオイル及びガスの位置を識別するので非常に重要である。残りのオイル又はガスの分布に関する知見は、さらなるフラッディング/処理の候補と考えられ得るのか否かを特定する上で助けとなる。その理由は、残りのオイル及びガスは、多数の既知の候補となる手法（たとえば第三次フラッディング、低塩分フラッディング、化学フラッディング、界面活性剤の注入、ガスの注入等）によって移動可能であるためである。これらの手法による流体の移動もまた、この方法によって直接的に検討することができる。
【０１１７】
図32Aは3D断層画像から得られた同一の薄片を図示している。前記3D断層画像は、様々な濡れ性条件下フラッディングを行った結果得られた様々な残留炭化水素分布を表している。3210は、水が濡れる条件下で水：オイル系へ水をフラッディングさせた後の残留炭化水素分布を示している。3212は、混合物の濡れ性(MW)条件下で水：オイル系へ水をフラッディングさせた後の残留炭化水素分布を示している。流体の分布はそれぞれ異なり、かつ定量化可能である。たとえば図32Bのグラフ3214乃至3220は、これらの様々な濡れ性条件下での残留炭化水素の小塊のサイズ（水のフラッディング後でのオイルを標的とする上で重要である）変化を定量化することができる。図32Bでは、残留炭化水素の小塊のサイズは、自発的吸水(SI)及び強制吸水(FI)の後に、水が濡れる条件（グラフ3218）から混合物が濡れる条件（グラフ3220）で変化する。
【０１１８】
図33は、どのようにして3D画像レジストレーションが、多数の長さスケールでの試料の特性を特定するのに利用可能であるのかを図示している。多くの材料は、（nmからcmまで）様々な寸法の3D構造部位を含んでいる。一の単一3D画像では、複数スケールの全てで全情報を捉えることができない（たとえばそのような画像は100000003すなわち200兆バイト（200テラバイト）を必要とする）。多数の長さスケールを探索するのに必要なギャップの架橋を助けるため、撮像は複数のスケールで行われて良く、かつ、定量的なレジストレーションは、大きな領域上での微細構造のマッピングに用いられて良い。
【０１１９】
mmからnmスケールにまでわたる部位を有する試料の画像の例が図33に図示されている。薄片3310によって表される3D画像が、20μmのボクセルサイズにて4cm直径の視野で得られる。これにより、20μmの材料部位と4cmの材料部位との識別が可能となる。しかし高解像度で部位を探索するため、4cmのスケールで本来撮像された試料の1つ以上の部分集合が再撮像されなければならない。この画像は本来の画像3310の部分集合で、かつ2.5μmの解像度で得られた。ここで部位は、改善された解像度で識別可能である。前記改善された解像度は元の解像度の約8倍よりも良好である。これら2つの画像から得られた情報は、3Dの画像データの組を得ることによって合成されて良い。前記3Dの画像データの組は、画像3310の各ボクセルに対してマッピングされた画像3312の高解像度画像データを有する。
【０１２０】
図34は、大きな画像3410へ埋め込まれる高解像度画像3412を図示している。高解像度画像3412は、画像3410内での部位の識別を容易にする。さもなければこれは、多数の複製後での「解像」が難しくなると思われる。
【０１２１】
図35は、各異なる解像度を有する複数の画像の同一の再分割した領域を並べて比較したものを図示している。低解像度の画像3510では解像されない部位が高解像度画像3512では識別することができる。画像レジストレーションの厳密かつ定量的な性質を有するため、高解像度の画像で識別された部位を利用することで、大きな画像中の不十分すなわち解像できない部位を「加える」ことができる。画像2510中で識別されるグレースケールは、高解像度画像3512内で明確に識別される相へ直接マッピングされて良い。たとえば画像3512内で明確に観察される有孔率と孔のサイズは、左側の画像中で識別されるグレースケール値に相関付けられて良い。図23及び図24に示された3D体積への2Dデータの伝播/合成と同様に、3412と3512から得られる関連する3D構造の情報（たとえば孔のサイズ、孔の形状、グレイン構造）は、より大きな体積（それそれ3410と3510）で画像データ全体にわたって伝播されて良い。
【０１２２】
F. 他の態様
開示された方法及びシステムは、3DのX線マイクロCTデータへの、（鉱物学及び探索のサブミクロンスケールである）従来の顕微鏡手法から取得可能な詳細な構造の情報の結合を助ける。特に、開示された方法は、マイクロCTから得られた3D画像への詳細な2D顕微鏡画像の最適なレジストレーションを可能にする。この結果、従来のマイクロCTイメージングの能力は顕著に拡張する。
【０１２３】
それに加えて、様々な濡れ性条件及び飽和状態下における同一の3D試料上での3D孔スケール分布を測定する方法並びにシステムが開示されている。当該方法は、ビームラインから有孔性試料を除去する手順、実験室内で系についての一組の実験を行う手順、及び、その後前記試料を再撮像することで、3D画像のボクセル−ボクセルのレジストレーションを行う手順を可能にする。
【０１２４】
画像データを処理する開示された方法を適用する結果、試料の1つ以上の特性に関する情報を有する合成画像データが得られる。得られた情報は、各独立した画像から得ることのできる情報に追加される。さらに情報は画像間の重なり領域を超えて補間されて良いので、当該方法は、レジストレーション処理に関与する本来の画像から得ることのできる結合された情報に追加される情報を得ることを可能にする。
【０１２５】
さらに当該方法は、試料の撮像領域全体にわたる−ひいては試料全体にわたる−その試料の特性の空間分布を推定することを可能にする。ここで「特性」という語は、非常に広い意味で用いられていて、かつ試料の情報を有する如何なる特徴をも含むと解される。そのような特性の例は以下の1つ以上を含んで良い。それは、物理的特性、構造特性、化学特性、鉱物学的特性、透水率、音響特性、グレインサイズ、グレイン境界、セメント相、地質学上外見の種類と地質学上の岩石の種類、鉱物相の構造、又は試料の孔への流体分布である。特に開示された方法のこのような有利な能力は、2つの3D画像の処理のレジストレーションの場合にも適用可能である。
【技術分野】
【０００１】
撮像手法、特に画像データを処理することで、その画像データが取得される物理的試料の特性の知見を得る方法及びシステムに関する。
【背景技術】
【０００２】
地球科学におけるコア材料の画像解析は、基本的には堆積学者及び記載岩石学の分野であった。光学顕微鏡及び電子顕微鏡(SEM)から得られる2次元(2D)画像データは、多孔率、孔のサイズ、グレインサイズ、フローユニット、透水率、速度、及び圧縮率の推定に利用されてきた。従来の記載岩石学の手法によって、鉱物の相情報及び鉱物の相の起源（たとえば摩滅又は自生）を特定することが可能となった。ナノメートル長さのスケールの探索が可能である。
【０００３】
他方、ミクロンスケールの孔構造の詳細な3次元(3D)画像を生成する能力を有するX線断層診断が近年、十分に確立された2D顕微手法にとって有用かつ不可欠なものとして受け入れられてきた。高品質の「ターンキー」X線断層診断システムが利用可能となったことは近年、これらのシステム利用の急速な増大を後押ししてきている。これらのシステムによって、3次元における有孔性材料についての孔/グレイン情報を得ることが可能となる。残念なことに、従来のマイクロCTイメージングでは鉱物の識別が不十分で、従来のマイクロCTイメージングの空間解像度は約1μmに制限される。
【０００４】
コア材料の孔/グレイン/粘土スケールでの特性の評価及び測定への注目が近年高まっている。このスケールでの特性の理解は、石油及びガス産業における応用にとって重要である。コア試料の解析は、カギとなる油層岩石物理学上のフロー特性及び多相フロー特性を生成するのに用いられる。これらの特性は、石油会社が、油田及びガス田の発見、油田及びガス田を生産状態にすること、並びに油田及びガス田の操業の際に直面する高い財務上のリスクを緩和する上で重要である。新たな油田及びガス田を生産状態にすることのコスト又は既存の油田及びガス田の寿命を延ばすことから得られると思われる利益を比較すると、解析それ自体のコストは低い。コアの解析は相反するデータ並びに解釈及び再現の困難なデータをよく与える恐れがあるという事実にもかかわらず、コアの解析は、貯蔵量の推定及び回復速度の予測をする上で業界標準データとなる。そのような困難には、多相フロー特性をよりよく理解するため、基本レベルで解明される必要がある複雑な界面現象が少なくとも部分的に起因している。
【０００５】
コアの実験室内での従来の多相フロー実験において得られた測定結果は、岩石の孔スケール構造と、及び、流体/流体と流体/岩石の相互作用の界面特性の両方の研究に用いられている。多相フローへ適用される孔レベルでの変位のモデル化規則及び方法の開発に大きな関心が存在する。様々な濡れ性のシナリオの下（水で濡れている、混合物で濡れている、大きい/小さい）での、個々の孔内での流体占有についての規則を含む孔ネットワークモデルが、現実の有孔性材料における多相フロー特性の理解を改善させる試みとして発展している（たとえば非特許文献1及び非特許文献2を参照のこと）。今日、様々な流体相の孔ネットワークモデルによる記述を、現実の濡れ性条件下での流体相の分布に関する実験によって得られた孔レベルの情報に対して直接校正するのは不可能である。従って孔のスケールでのモデルを直接孔のレベルで校正するのは不可能である。
【０００６】
CT画像に関する直接的なシミュレーションが、単相の有孔性材料の特性−たとえば透水率、伝導率、及び水銀注入時の毛細管圧力曲線−を予測するのに利用可能であることが示された（たとえば非特許文献3を参照のこと）。
【０００７】
これまでの研究は、マイクロCT画像化実験に基づいて3Dでの孔スケールでの流体分布を特定することが可能であることを示してきた（たとえば非特許文献4を参照のこと）。しかしこれらの研究は、X線CTのビームラインから試料を取り除くことなく、検査中にコア試料をフラッディングさせる必要があることに深刻に悩まされている。このため、3Dにおいて孔スケールで行われる実験は、単純なフラッディング実験に制限されている。
【０００８】
しかもこれらの実験は、長い取得期間にわたってコア材料の位置をミクロンスケールで完全に維持する必要があることで、制約を受けてきた。よって、平衡に達するのにかなりの回数を含む実験（たとえば、油層条件下での定常状態の相対的な透水性測定、天然の原油及びブラインにおけるコアの劣化、有孔性板のフラッディング等）が、数週間またさらには数ヶ月、試料をCT装置で維持するのに必要となるだろう。
【先行技術文献】
【非特許文献】
【０００９】
【非特許文献１】Morrow & Mason, "Recovery of Oil by spontaneous imbibition", Current Opinion in Colloid & Interface Science, Vol. 6, pp. 321-337(2001)
【非特許文献２】H. Behbahani and M. Blunt, "Analysis of Imbibition in mixed wet rocks using pore scale modelling" , SPE 90132, presented at SPE Annual Technical Conference, Houston, 2004
【非特許文献３】C. H. Arns, F. Bauget, A. Ghous, A. Sakellariou, T. J. Senden, A. P. Sheppard, R. M. Sok, W. V. Pinczewski, J. Kelly, and M. A. Knackstedt, "Digital core laboratory: Petrophysical analysis from 3D imaging of reservoir core fragments ", Petrophysics, 46(4), 260-277, 2005.
【非特許文献４】Seright et al., "Characterizing disproportionate permeability reduction using synchrotron X-ray computed tomography", SPE Reservoir Evaluation and Engineering, October 2002, pp. 355-364
【非特許文献５】C. H. Arns, A.P. Sheppard, R. M. Sok and M. A. Knackstedt, NMR petrophysical predictions on digitized core images Petrophysics, 48 (3), 202-221 (2007)
【非特許文献６】M. Knackstedt, C. Arns, A. Sheppard, T. Senden, R. M. Sok, W. Pinczewski and M. Archie 's exponents in complex lithologies derived from, 3D digital core analysis, Ioannidis, Society of Petrophysicists and Well Log Analysts 48th Annual Logging Symposium, Paper Z, June 4-6, 2007, Austin, USA.
【非特許文献７】Morrow, "Wettability and its effect on oil recovery", SPE 21621, published in Journal of Petroleum Technology, December 1990, p.1476
【非特許文献８】Graue et al., "Alteration of wettability and wettability heterogeneity" Journal of Petroleum Science & Engineering, Vol. 33, pp. 3-17
【非特許文献９】M. Prodanovic et al "3D image-based characterization of fluid displacement in a Berea core" Advances in Water Resources 30, 2007, p214
【非特許文献１０】Dautriet et al "Laboratory determination of stress path dependency of directional permeability of Estaillades limestone", Symposium of the Society of Core Analysts, Abu Dhabi, November 2008, paper SCA2008-26
【発明の概要】
【発明が解決しようとする課題】
【００１０】
本発明の目的は、既知の装置の1つ以上の欠点を実質的に克服する、若しくは少なくとも改良すること、又は、有用な代替解決策を提供することである。
【課題を解決するための手段】
【００１１】
本開示の第1態様によると、試料の画像データを処理する方法が供される。当該方法は、前記試料の少なくとも一部が重なっている複数の空間領域の第1画像と第2画像を登録する手順、及び、前記の登録された画像からのデータを処理することで、前記第1像と第2像から得られる画像に加えられる前記試料についての情報を有する合成画像データを取得する手順を有する。
【００１２】
本開示の第2態様によると、前記第1態様の方法に係る工程を実施することが可能なコンピュータプログラムが供される。
【００１３】
本開示の第3態様によると、前記第2態様のコンピュータプログラムを有するコンピュータによる読み取り可能な媒体を有するコンピュータプログラム製品が供される。
【００１４】
本開示のさらなる態様によると、試料の画像データを処理するエレクトロニクスシステムが供される。当該システムは、前記試料の少なくとも部分的に重なっている複数の空間領域の第1画像と第2画像を取得するように備えられた撮像手段、前記の取得された画像を記憶するように備えられた記憶装置、並びに、前記の取得された画像データを取得及び登録することで、前記試料についての情報を有する合成画像データを得る工程であって、前記情報は前記第1画像と第2画像から得られる像に加えられる、工程、を有する。
【００１５】
本開示の他の態様もまた開示されている。
【図面の簡単な説明】
【００１６】
【図１】本発明の一の実施例による方法の概略的なフローダイアグラムを図示している。当該方法は、2D及び3D画像を登録し、その後両画像からの情報を合成することに関する。
【図２Ａ】2D電子顕微鏡像を図示している。
【図２Ｂ】3D断層画像を図示している。
【図３】図2A及び図2Bから得られた登録された画像のアニメーションの最終フレームを図示している。
【図４Ａ】2D光学顕微鏡像を図示している。
【図４Ｂ】3D断層画像を図示している。
【図５】2D顕微鏡像を図示している。
【図６】登録された3D断層画像を図示している。図6の画像に係る解像度は図5の画像に係る解像度とは異なる。
【図７】2D画像から取得された孔構造の情報を3D画像データに合成することによって得られた合成画像を図示している。
【図８】ワーピング変換（warp transformation）が用いられた規則的な正方形グリッドの走査型電子顕微鏡(SEM)像を図示している。
【図９】ワープ変換ステンシルの概略図を示している。
【図１０】ワープ変換において画像の中心に最も近い5つの最大相関値を図示している。
【図１１Ａ】登録前に逆ワープされなかった画像から得られた像を図示している。
【図１１Ｂ】登録前に逆ワープされた画像から得られた像を図示している。
【図１２】本発明の第2実施例による方法の概略的フローダイアグラムを図示している。当該方法は、3D画像同士の位置合わせを行い、その後両画像からの情報を合成することに関する。
【図１３】流体吸入前の試料を図示している。
【図１４】流体吸入後の試料を図示している。
【図１５】試料内部での流体分布を示す動画からのスナップショットである。前記流体分布は3D画像を位置合わせすることによって得られた。
【図１６Ａ】3D相分布を図示している。
【図１６Ｂ】孔のサイズ分布と流体含有量との相関を示している。
【図１７】本発明の第1実施例に係る計算手法のアルゴリズムを図示している。
【図１８】本発明の第2実施例に係る計算手法のアルゴリズムを図示している。
【図１９】図17と図18に示されたアルゴリズムの手順2のフローダイアグラムの詳細を図示している。
【図２０】3D断層画像内部に登録された2D顕微鏡像を図示している。
【図２１】2つの登録された画像中の部位が等しいことを示している。前記2つの登録された画像のうちの一は2DのSEM像で、前記2つの登録された画像のうちの他は3D断層画像から得られた対応部分の薄片である。
【図２２】2DのSEM像と、該2DのSEM像と等価な、登録された3D断層画像から得られた薄片を図示している。
【図２３】3D断層画像から得られた減衰データと2D顕微鏡像から得られた多孔率との間での依存性を表す関数のプロットである。
【図２４Ａ】解像度の異なる2つの画像と、該2つの画像から得られる合成データに基づく画像を図示している。
【図２４Ｂ】図24Aの画像から得られる画像データに基づく、局所的な孔のサイズ分布（Y軸）とX線減衰から得られたミクロ孔マッピングとの相関を示している。
【図２５】3D断層画像から得られた減衰データに対する、2Dミクロ孔から得られた鉱物相の情報のマッピングを図示している。
【図２６】試料から3D画像を取得する手順、前記試料を撮像装置から取り外して、前記試料の実験を行う手順、前記試料を再撮像する手順、及び、前記実験の前後で得られた2つの3D画像を登録することで、前記試料における構造変化又は他の変化を特定する手順の概念を図示している。
【図２７】試料から3D画像を取得する手順、前記試料を撮像装置から取り外して、前記試料の実験を行う手順、前記試料を再撮像する手順、及び、前記実験の前後で得られた2つの3D画像を登録することで、前記試料における構造変化又は他の変化を特定する手順の概念を図示している。
【図２８】カルボン酸で処理された試料の登録3D画像の2つの対応する薄片を図示している。明らかな差異は、試料材料特性への酸の効果に起因する。
【図２９】試料の一連の状態から決まる一連の曲線を図示している。各曲線は、特定の値の水飽和圧力及び毛管圧力によって特徴付けられる。
【図３０】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３１】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３２Ａ】図29のグラフに関する、様々な状態をとる試料の登録された画像を図示している。
【図３２Ｂ】様々な濡れ性条件下での残留炭化水素の小塊のサイズ変化を定量化している。
【図３３】複数の解像度を有する3D画像を登録する応用例を図示している。
【図３４】高解像度の3D画像を低解像度の3D画像へ埋め込む応用例を図示している。
【図３５】直接比較するために配置された異なる解像度の登録3D画像の対応する薄片を図示している。
【図３６Ａ】記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【図３６Ｂ】記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【発明を実施するための形態】
【００１７】
ここで1つ以上の本発明の実施例について図面を参照しながら説明する。
【００１８】
図1乃至図35は、登録された2D画像と3D画像を合成することで得られたデータを処理するための開示された方法に係る様々な態様、及び当該処理方法によって可能となる実際の用途を示している。
【００１９】
A ハードウエアの実装
図36Aと図36Bは一団をなすように、記載された構成の実施が可能な汎用コンピュータシステムの概略的ブロック図を形成する。
【００２０】
図36Aから分かるように、コンピュータシステム3600は、コンピュータモジュール3601、入力装置−たとえばキーボード3602、マウスポインタ装置3603、スキャナ3626、及びマイクロホン3680−、及び、出力装置−プリンタ3615、ディスプレイ装置3614、及びラウドスピーカー3617を含む−によって形成される。1つ以上の撮像装置3627もまたコンピュータモジュール3601と直接接続して良い。これらの装置は、後述の方法によって処理される画像データを供する。
【００２１】
外部の変調器−復調器（モデム）トランシーバ装置3616は、接続3621を経由してコミュニケーションネットワーク3620とのやり取りを行うため、コンピュータモジュール3601によって使用されて良い。ネットワーク3620は広域ネットワーク(WAN)−たとえばインターネット又はプライベートWAN−であって良い。接続3621が電話線である場合、モデム3616は従来の「ダイアルアップ」モデムであって良い。あるいはその代わりに、接続3621が高容量（たとえばケーブル）接続である場合、モデム3616はブロードバンドモデムであって良い。ネットワーク3620へのワイヤレス接続用にワイヤレスモデムが用いられても良い。
【００２２】
コンピュータモデム3601は一般的に、少なくとも1つのプロセッサ、及び、たとえば半導体ランダムアクセスメモリ(RAM)及び半導体読み出し専用メモリ(ROM)から構成されるメモリユニット3606を有する。図36Aに図示されたコンピュータシステムは2つのプロセッサ3605Aと3605Bを有する。2つのプロセッサ3605Aと3605Bは、後述する方法の手順の少なくとも一部において実装される並列処理を可能にする。モジュール3601はまた多数の入出力(I/O)インターフェースをも有する。前記多数の入出力(I/O)インターフェースは、ビデオディスプレイ3614、ラウドスピーカー3617、及びマイクロホン3680と結合するオーディオ−ビデオインターフェース3607、キーボード3602、マウス3603、スキャナ3626、撮像装置3627、及び任意のジョイスティック（図示されていない）用のI/Oインターフェース3613、並びに、外部モデム3616及びプリンタ3615用のインターフェース3608を有する。一部の実施例では、モデム3616はコンピュータモジュール3601内部−たとえばインターフェース3608内部−に組み込まれて良い。コンピュータモジュール3601はまたローカルネットインターフェース3611をも有する。ローカルネットインターフェース3611は、接続3623を介して、コンピュータシステム3600を、ローカルエリアネットワークとして知られているローカルコンピュータネットワーク3622へ接続することを可能にする。また図示されているように、ローカルネットワーク3622はまた、接続3624を介して広域ネットワーク3620へも接続して良い。広域ネットワーク3620は一般的には、所謂「ファイアーウォール」デバイス又は類似の機能を有するデバイスを含む。インターフェース3611は、Ethernet（商標）回路カード、Bluetooth（商標）ワイヤレス構成、又はIEEE802.11ワイヤレス構成によって形成されて良い。
【００２３】
インターフェース3608と3613はシリアル接続及び/又はパラレル接続が可能である。シリアル接続は一般的に、ユニバーサルシリアルバス(USB)標準規格に従って実装されて良く、かつ対応するUSBコネクタ（図示されていない）を有する。他の記憶装置−たとえばフロッピーディスクドライブ及び磁気テープ装置（図示されていない）−が用いられても良い。光学ディスクドライブ3612は一般的には、不揮発性データ源として機能するように供される。持ち運び可能なメモリ装置−たとえば光学ディスク（たとえばCD-ROM,DVD）、USB-RAM、及びフロッピーディスク−が、システム3600への適切なデータ源として用いられて良い。
【００２４】
コンピュータモジュール3601の構成部品3605乃至3613は一般的に、相互接続バス3604を介して、当業者に既知のコンピュータシステム3600が従来の動作モードとなるようにやり取りを行う。上述の構成を実施することが可能なコンピュータの例には、IBM-PC、サンマイクロシステム(Sun)のスパークステーション(Sparkstations)、アップル(Apple)社のマッキントッシュ(Macintosh)（商標）又はこれらの進化型が含まれる。
【００２５】
上述の画像処理方法はコンピュータシステム3600を用いて実装されて良い。これから説明する図1乃至12の手順及び図17乃至19の計算アルゴリズムは、コンピュータシステム3600内部で実行可能な1つ以上のソフトウエアアプリケーションプログラム3633として実装されて良い。特に、上述の方法の手順は、コンピュータシステム3600内部で実行されるソフトウエア3633内の命令によって行われる。ソフトウエアの命令3631は1つ以上のコードモジュールとして生成されて良い。各コードモジュールは特定のタスクを実行することを意図している。そのソフトウエアはまた2つの別個の部分に分割されて良い。第1部分と対応するコードモジュールは当該方法の手順を実行する。第2部分と対応するコードモジュールは前記第1部分とユーザーとの間のユーザーインターフェースを管理する。
【００２６】
ソフトウエア3633は一般的に、コンピュータによる読み取り可能な媒体からコンピュータシステム3600へロードされ、続いて典型的には、図36Aに図示されているようにHDD3610に記憶されるか、又はメモリ3606に記憶される。その後ソフトウエア3633はコンピュータシステム3600によって実行可能となる。一部の例においては、アプリケーションプログラム3633は、ユーザーへ供給され、1つ以上のCD-ROM3625上でエンコードされ、かつメモリ3610又は3606に記憶される前に、対応するドライブ3612を介して読み出されて良い。あるいはその代わりにソフトウエア3633は、ネットワーク3620又は3622からコンピュータシステム3600によって読み出されて良いし、又は他のコンピュータによる読み取り可能な媒体からコンピュータシステム3600へロードされても良い。コンピュータによる読み取り可能な媒体とは、実行及び/又は処理のためコンピュータシステム3600へ命令及び/又はデータを供することに関与する如何なる記憶媒体をも指称する。そのような記憶媒体の例には、フロッピーディスク、磁気テープ、CD-ROM、ハードディスクドライブ、ROM、集積回路、USBメモリ、磁気光学ディスク、又はPCMCIAカードのようなコンピュータによる読み取り可能な媒体等が含まれる。これらはコンピュータモジュール3601の内部又は外部であるかを問わない。コンピュータモジュール3601へのソフトウエア、アプリケーションプログラム、命令、及び/又はデータの提供にも関与することが可能なコンピュータによる伝送可能な媒体の例には、高周波又は赤外伝送チャネル及び他のコンピュータ若しくはネットワーク装置へのネットワーク接続、並びに、eメール伝送及びウエブサイトに記録された情報を含むインターネット及びイントラネット等が含まれる。
【００２７】
上述のアプリケーションプログラム3633の第2部分及び対応するコードモジュールは、ディスプレイ3614上に描画又は表示される1つ以上のグラフィカルユーザーインターフェース(GUI)を実装するように実行されて良い。典型的にはキーボード3602とマウス3603の操作を介して、コンピュータシステム3600とアプリケーションのユーザーは、機能が一致するように前記（複数の）GUIに係るアプリケーションへの制御コマンド及び/又は入力を供するようにしてインターフェースを制御することができる。機能が一致するような他の形式のユーザーインターフェースが実装されても良い。そのようなインターフェースとはたとえば、ラウドスピーカー3617を介したスピーチプロンプト出力を利用するオーディオインターフェース及びマイクロホン3680を介したユーザー音声コマンド入力である。
【００２８】
図36Bは、プロセッサ3605A又は3605Bのうちのいずれか1つ及び「メモリ」3634の詳細な概略ブロック図である。メモリ3634は、図36Aのコンピュータモジュール3601によってアクセス可能な全メモリ装置の論理的集合を表す（HDD3610及び半導体メモリ3606も含まれる）。
【００２９】
コンピュータモジュール3601は最初電源がオンにされるとき、ハードウエアの初期化と動作チェック(POST)プログラム3650が実行される。POSTプログラム3650は一般的には半導体メモリ3606のROM3649内に記憶される。ハードウエア装置−たとえばROM3649−のような恒久的に記憶されるプログラムはファームウエアと呼ばれることがある。POSTプログラム3650は、適切に機能することを保証するためにコンピュータモジュール3601内部でハードウエアを点検し、かつ典型的には、プロセッサ3605A又は3605Bのうちの対応する1つ、メモリ3609又は3606、及び、校正操作を行うために一般的にはROM3649内に記憶される基本入出力システムソフトウエア(BIOS)モジュール3651をチェックする。一旦POSTプログラム3650が問題なく動作すると、BIOS3651はハードディスクドライブ3610を起動する。ハードディスクドライブ3610が起動することで、該ハードディスクドライブ3610内に存在するブートストラップローダープログラム3652が、プロセッサ3605A又は3605Bのうちの対応する1つによって実行される。これにより、OS3653はRAMメモリ3606へロードされ、RAMメモリ3606上でOS3653は動作を開始する。OS3653は、様々な高次の機能を実現するための、プロセッサ3605A又は3605Bのうちの対応する1つによって実行可能なシステムレベルのアプリケーションである。様々な高次の機能には、プロセッサ管理、メモリ管理、デバイス管理、記憶装置管理、ソフトウエアアプリケーション管理、及び汎用ユーザーインターフェース管理が含まれる。
【００３０】
OS3653は、コンピュータモジュール3601上で動作するプロセス又はアプリケーションの各々が十分なメモリを有することを保証するため、他のプロセスに割り当てられたメモリと相反することなく実行できるようにメモリ3609又は3606を管理する。さらに各プロセスが効率的に動作できるように、システム3600において利用可能な様々な型のメモリが適切に用いられなければならない。従って、メモリの一団3634は、どのようにしてメモリの特別な部分が割り当てられるのか（指定されない場合）を示すことを意図しておらず、むしろコンピュータシステム3600によってアクセス可能なメモリの一般的な視点及びそのようなメモリがどのようにして用いられるのかを示すことを意図している。
【００３１】
プロセッサ3605Aと3605Bの各々は多数の機能モジュールを有する。前記多数の機能モジュールには、制御ユニット3639、算術論理ユニット(ALU)3640、及び局所又は内部メモリ3648−キャッシュメモリと呼ばれることもある−が含まれる。キャッシュメモリ3648は一般的に、レジスタ部内に多数の記憶レジスタ3644-3646を有する。1つ以上の内部バス3641はこれらの機能的モジュールと機能的に相互接続する。プロセッサ3605Aと3605Bの各々は一般的に、接続3618を用いて、システムバス3604を経由して外部装置とのやり取りを行う1つ以上のインターフェース3642をも有する。
【００３２】
アプリケーションプログラム3633は一連の命令3631を有する。一連の命令3631は、条件的分岐及び巡回命令を有して良い。プログラム3633はまた、プログラム3633の実行に用いられるデータをも有して良い。命令3631及びデータ3632は、メモリ位置3628-3630及び3635-3637にそれぞれ記憶される。命令3631とメモリ位置3628-3630の相対的なサイズに依存して、メモリ位置3630内に図示された命令によって示されているように、特定の命令が単一のメモリ位置に記憶されて良い。あるいはその代わりに、命令は多数の部分に分割されて良い。メモリ位置3628-3629内に図示されている命令セグメントによって示されているように、各部分は各独立したメモリ位置に記憶される。
【００３３】
一般的には、プロセッサ3605Aと3605Bには、内部で実行される1組の命令が与えられる。続いてプロセッサ3605Aと3605Bは、他の組の命令の実行によって応答するまで、後続の入力を待つ。各入力は多数の入力源のうちの1つ以上から供されて良い。前記入力には、1つ以上の入力装置3602と3603によって生成されるデータ、ネットワーク3620と3622のうちの1つを経由して外部から受信されるデータ、記憶装置3606と3609のうちの1つから取得されるデータ、又は、対応する読み取り装置3612へ挿入された記憶媒体3625から取得されるデータが含まれる。1組の命令が実行される結果、場合によっては、データが出力されることがある。実行には、メモリへのデータ又は変数の記憶もが含まれて良い。
【００３４】
本明細書において以降に開示される計算機の構成は、メモリ3634内の対応するメモリ位置3655-3658内に記憶される入力変数3654を用いる。実行された計算は、メモリ3634内の対応するメモリ位置3662-3665内に記憶される出力変数3661を生成する。中間変数が、メモリ位置3659、3660、3666、及び3667内に記憶されても良い。
【００３５】
動作するプロセッサ3605A並びに/又はプロセッサ3605Bのレジスタ部3644-3646、算術論理ユニット(ALU)3640、及び制御ユニット3639は、協働することで、プログラム3633を構成する命令の組に含まれる命令毎に「取得、復号、及び実行」のサイクルを行う必要があるミクロ動作のシーケンスを行う。各取得、復号、及び実行サイクルは、
(a) メモリ位置3628からの命令3631の取得又は読み取りを行う取得動作、
(b) 制御ユニットがどの命令が取得されたのかを判断する復号動作、
(c) 制御ユニット3639及び/又はALU3640が前記命令を実行する実行動作、
を有する。
【００３６】
その後、次の命令のためにさらなる取得、復号、及び実行サイクルが実行されて良い。同様に、制御ユニット3639がメモリ位置3632へ値を記憶又は書き込みするための記憶サイクルも実行されて良い。
【００３７】
その後、図17-19のアルゴリズムに含まれる各手順又は下位の処理が、プログラム3633の1つ以上の部分に関連付けられ、かつ、プロセッサ3605A（及び/又は3605B）内のレジスタ部3644-3647、ALU3640、及び制御ユニット3639によって実行される。プロセッサ3605A（及び/又は3605B）内のレジスタ部3644-3647、ALU3640、及び制御ユニット3639は、協働することで、プログラム3633を構成する命令の組に含まれる命令毎に「取得、復号、及び実行」のサイクルを行う必要があるミクロ動作のシーケンスを行う。
【００３８】
あるいはその代わりに、画像処理の開示された方法は、様々な計算ルーチンの機能又は下位の機能を実行する専用ハードウエア−たとえば1つ以上の集積回路−で実装されても良い。そのような専用ハードウエアには、グラフィックプロセッサ、デジタル信号プロセッサ、又は1つ以上のマイクロプロセッサと付属のメモリが含まれて良い。
【００３９】
2つの画像のデータの組が有する情報の比較及び合成を可能にするため、各対応する画像データの組のうちの少なくとも1つは、2つの組が1つの座標系を表すように変換されなければならない。光イメージングの分野では、この変換は「レジストレーション（画像の位置合わせ）」と呼ばれ、かつ、1つの座標系を表すように「位置合わせ」される画像は「レジストレーションされた」画像と呼ばれる。
【００４０】
B レジストレーションされた2D画像と3D画像の処理
以降の開示は、3D撮像−たとえば従来のマイクロCT−の能力を拡張し、かつ2D顕微手法の範囲で得られる相補的データ及び/又は高解像度データを用いることによる、3D画像データの校正を可能にする方法に関する。そのようにして得られた合成データの一の用途は、コア材料/試料の地質学的評価である。具体的には、開示された方法は、3D画像上及び該3D画像内部で2D画像の画像データのレジストレーションに関する。
【００４１】
3D画像体積内部に存在する2D顕微鏡画像のピクセル−ボクセル識別を3D画像と同程度の解像度で行うことで、3D画像中のX線グレースケールレベルを顕微手法から得られたある範囲の特性にマッピングすることが可能となる。また前記ピクセル−ボクセル識別は、鉱物学、元素組成、表面エネルギー、及び音響/弾性特性のうちの少なくとも1つを含む特性の研究に用いられる。従来のX線マイクロCTが約1μmの解像度を有する一方で、顕微手法は、部位及び特性の探索をナノスケールで行うことを可能にする。試料内部の少なくとも一部が重なっている複数の空間領域の3D画像に対して2D画像を位置合わせして、該2D画像を3D画像へマッピングする能力は、現時点においてX線CT解析では微細すぎると考えられている特性についての情報の取得を可能にする。レジストレーションされ、かつ処理された画像は、任意に配向し、スケールを変換し、かつ相互に並進して良い。さらに2つの3D画像の処理のレジストレーションが行われる場合−この場合2つのレジストレーションされた3D画像が任意に配向し、スケールを変換し、かつ相互に並進して良い−にも、開示されたこの有利な能力はまた適用可能である。2D画像と3D画像との間でのレジストレーションの場合では、部分的な重なりは、3Dイメージング手法によってイメージングされる3D体積と、2Dイメージングによって試料上又は該試料内部にイメージングされる任意に配向した面との間での線又は面の重なりであることに留意して欲しい。3D-3Dレジストレーションの場合では、部分的な重なりは一般的には、1つ以上の各対応するイメージング手法によってイメージングされる、試料の2つ以上の3D領域間での空間的な重なりである。
【００４２】
開示された方法によって実行されるマッピングは、ある範囲の2D顕微手法によって得られた従来の堆積データ（たとえばコアの記載、記載岩石学、地質学上の外見、地質学上の岩石の種類、鉱物の相、鉱物の起源、弾性/音響特性、表面化学に関連するデータ）を、マイクロCT解析によって得られる3D画像データ（たとえば減衰、グレイン形状、グレインの接続性、孔の接続性、孔の形状に関連するデータ）に相関させることを可能にする。顕微手法では、ナノメートルでの探索が可能である。これはマイクロCTの解像度よりも小さい。続いて得られた高解像度データは、3DマイクロCT画像からのX線減衰データに直接相関付けられて良い。よって従来の堆積学的情報は、油層岩石物理学的特性（透水性、排水毛管圧力、音響及び地震特性、抵抗力、NMR）及び油層の設計に係る特性（たとえば相対的透水性、残留飽和度、吸水毛管圧力）の予測に関連するデータに相関付けることができる。
【００４３】
3D画像に対する2D画像のレジストレーションから得られる合成情報
当該方法の詳細が図1に図示されている。一旦有孔性材料が取り出されると、手順102では、マイクロ断層撮像によって3D画像が得られる。試料から得られた1つ以上の2D部分又は露出した平面は、手順104において、顕微鏡観察/分光用に収集及び調製される。コア材料の残りは品質制御実験（たとえばヘリウム及び水銀注入による多孔度測定）に用いられて良い。
【００４４】
顕微鏡観察/分光は、手順106において、以下の手法のうちの1つ以上によって行われる。
− 走査電子顕微鏡(SEM)：2次電子及び後方散乱電子
− 光学顕微鏡
− 走査音波顕微鏡
− レーザー共焦点顕微鏡
− 集束イオンビーム−走査電子顕微鏡(FIBSEM)
− エネルギー分散分光(EDAX)
− 2次イオン質量分光(SIMS)
− 分光学的手法：赤外、UV−可視、及びラマン分光が含まれるが、これらに限定されるわけではない
− 他の顕微又は分光手法
得られた顕微/分光2D画像データは、従来手法で処理されることで、ある範囲の長さスケール（nm〜μm）にわたって、2Dにおいて試料の様々な構造の特性が取り出される（手順110）。この範囲であれば、たとえば鉱物学、表面特性、局所的な表面相互作用、音響特性等の特性の解析は可能である。例として、光学顕微鏡による標準的な記載岩石学的解析は、以下のうちの少なくとも1つに係る特性の測定を可能にする。
− 有孔性
− 鉱物（たとえば石英、長石、方解石、粘土、黄鉄鉱等）
− グレイン境界、セメント相の特定
− 地質学的外見の種類の特定
− 地質学上の岩石の種類
当該方法の手順108は、3D画像上の1つ以上の2D画像のピクセル−ボクセルレジストレーションを実行する。手順110は従来の2D画像解析（たとえば分光写真解析）を実行する。これにより、当該方法の手順112において、定量的な合成データを得ることが可能となる。手順112は、1種類の手法によって得られたデータに値を追加することが可能である。合成された情報の一の用途は、高解像度の高品質2D顕微鏡解析データとの比較によって、3D画像データ上での直接的な品質制御を行うことである。合成データはまた、高解像度2D顕微鏡画像から得られた孔及び鉱物相の情報を、サブCT解像度又は従来のCTの解像度を有する3D画像から得られたデータに直接相関付けるのに用いることもできる。そのような相関は以下のうちの1つに用いられて良い。
− マイクロ多孔性測定の直接検査
− サブCT解像度スケールから得られたグレイン/孔のサイズ情報のマイクロ多孔性データへの追加
− 2D顕微手法及び分光手法によって得られた情報（たとえば鉱物学、膠着(cementation)、続成作用の履歴、高解像度（マイクロ）多孔性、実効透水率、グレインサイズ、サブミクロン解像度での孔のサイズ、音響特性、及び表面化学特性）の、マイクロCTから得られた3D画像データへの追加
また上述のデータ補正を用いることによって、ある範囲の物理的特性を任意の手順114においてシミュレーションすることができる。そのような物理的特性には、様々な流体相の相対透水率、誘電応答、NMR応答、及び音響地震応答が含まれる。典型的な適用可能なシミュレーションについての情報は、特許文献3、特許文献5、及び特許文献6で見つけることができる。
【００４５】
当該手法の応用例
1. 一例では、3D画像データの直接的な品質制御は、該3D画像データを標準的な2D顕微鏡解析から得られたデータと比較することによって実現されて良い。図2及び図3は、岩石試料のピクセル−ボクセルレジストレーションを表している。前記岩石試料の孔のほとんどはマイクロCTデータの解像度で十分識別されうる。図2AはSEM画像を図示している。図2Bは、手順108（図1）に示されているように、開示された方法によってレジストレーションがなされた各対応する断層薄片を図示している。2D画像と3D画像との直接的な対応が存在する。如何なるばらつきも顕微鏡解析の試料調製に起因する。
【００４６】
2. 他の応用では、3D画像データからの減衰情報と、2D顕微鏡解析からの改善された情報とが直接的に相関付けられてよい。図4Aは、砂岩石試料の2D光学顕微鏡画像の一例を図示している。鉱物学上の広範囲にわたるばらつきを観測することができる。鉱物相分布がマッピングされて良い。ある範囲の古典的記載岩石学的解析が行われて良い。他方図4Bは、3D断層撮像から得られた対応するレジストレーションがなされた薄片を図示している。3D画像における減衰情報を、光学顕微鏡画像から得られた鉱物学上の情報に直接相関付けることに基づいて、顕微手法から得られた鉱物学上の情報は、3DマイクロCT画像データ内部の情報を追加/伝播するように相関付けられて良い。
【００４７】
別な例が図5に図示されている。図5は、50nm解像度で得られたミクロな試料の画像を図示している。図6は、2.5μmの解像度で得られた断層画像からのレジストレーションがなされた薄片を図示している。マイクロCTデータからのグレースケール減衰情報は、50nmスケール情報（たとえば高解像度での孔のサイズ、孔の形状、グレインのサイズ及び形状）に対して直接相関付けられて良い。
【００４８】
3. 当該方法はまた、様々な部位を識別する高解像度の顕微手法を用いる能力をも供する。前記様々な部位は、対応する3DマイクロCT画像と空間的に関連づけられる。一旦これらの部位が、空間的な関係によって、2D画像と3D画像との間での重なり領域について特定されると、残りの3D空間についてもこれらの部位を推定することができる。この手法は有効に顕微鏡データを利用することで、追加の情報を3D画像データに組み込む。図7A-7Cの画像は、どのようにして50nmスケールの2DのSEM画像から得られたグレインサイズと孔のサイズの情報が、図7Aの3D断層画像に埋め込まれた（マッピングされた）孔の構造の情報を生成することができるのかを示している。よって、図7AのマイクロCT画像からのデータと図7B及び図7Cの薄い部分のデータとを結合した合成画像データの組が生成される。そのような合成データは、前記薄い部分内部での微少孔領域の探索によって、ナノメートルサイズの孔のサイズ情報を抽出することを可能にする。nmスケール（50nm以下）で画像情報を探索することが可能であるため、これらのスケールでの3D画像における接続性の特性を測定することが可能となる。
【００４９】
一の応用では、孔のサイズとグレースケール減衰とが直接的に相関することで、本来は単純な有効率として割り当てられていた値を、ボクセルの孔のサイズ情報として3Dデータで伝播させることが可能となる。鉱物学上の情報、音響/弾性特性、及び表面化学情報もまた、最適となるように3DマイクロCTデータを追加する確率論的手法によってマッピングされて良い。
【００５０】
湾曲歪み補正：幾何学的画像補正
顕微鏡の中には、他の画像とのレジストレーションが行われる画像に幾何学収差を導入するものがある。このことはたとえば、画像中における2つの対象物間の測定間隔が、前記2つの対象物間の実際の間隔と一致しないことを意味すると考えられる。正確な定量的情報を推定する前、又は、これらの情報を適切に合成するために2つの画像が正確にレジストレーションされる前に、この幾何学な歪みは除去される必要がある。そのような歪みを画像から除去する処理は湾曲歪み補正と呼ばれる。
【００５１】
湾曲歪み補正は2つの段階からなる。最初に、観察方向に係るワーピング変換が決定される。続いて、逆ワーピング変換が画像に適用されることで、幾何学的に補正された画像が生成される。
【００５２】
ワーピング変換の一例は、正しい相対位置が既知である点からなる組で構成される画像を取得することである。よってこの情報は、歪んでいない画像中の点と歪んでいる画像中の点との間での滑らかなマップを構築するのに用いられる。任意の点の集合がこの目的に用いられて良い。しかし繰り返しになるが、補正の正確さは、その補正の複雑さと滑らかさに依存する。
【００５３】
以降では、2D画像のワーピング変換を自動で決定し得る典型的な方法について記載する。同じ原理は3Dの場合にも適用可能である。
【００５４】
湾曲歪み補正の第1段階は以下を含む。
1. 規則的な正方形グリッドの画像が取得される（図8）。この画像は、グリッド像と呼ばれ、画像中の歪みの決定に用いられる。
2. 前記グリッド像の中心に位置する小さな領域910−各方向において2,3のグリッドセルを含むのに十分な大きさで、かつ実線の暗い影付き四角形で表される−が抽出される。明るい影付き四角形912は、ステンシルパターンが繰り返される場所を示している。
3. 前記グリッド像と前記ステンシルとの間の相互相関が計算され、かつ相関ピークが決定される。前記相関ピークは、同一パターンが再度現れる場所を示している。小さな歪みについては、単純な並進で十分正確である。大きな歪みについては、前記ステンシルの回転及び伸縮が必要になる場合がある。いずれの場合でも、フーリエ変換が、これらの研鑽を効率的に行うのに用いられる。
4. 図10に図示されているように、前記像の中心に最も近い5つの相関ピーク並びに該ピークとグリッドの主軸u及びvとの関係が決定される。このことは、それ自体が中心位置であるだけではなく、前記グリッドの各主軸に沿った第1相関ピークでもある。これらの軸は前記画像の軸と位置合わせされている必要はない。
5. これら5つの相関ピークから、前記グリッドの主軸が決定されて良い。よってこれらの軸は、残りの相関ピークの正しい位置−つまり前記残りの相関ピークが歪んでいない画像中で存在したであろう位置−を推定するのに用いられる。
6. 正しい位置を各相関ピークに割り当てることで、この組についてのワーピング変換が定義される。よって、これらの点の間での位置についてのワーピング変換は適切な補間法によって定義される。
【００５５】
第2段階は前記画像の湾曲歪み補正を行う。この地点では、歪んでいない画像を得るため、前記歪んだ画像に逆ワーピング変換が適用される。つまり、前記歪んでいない画像中の各ピクセルは対応するワーピングされた位置へ移され、かつ前記ピクセルの値は前記歪んでいない画像中での補間によって求められる。前記歪んでいない画像に対して位置合わせされた規則的なグリッドによって処理することで、計算は大幅に単純化される。計算効率のため、この処理は2つのステップとして実行される。
1. 前記歪んでいない画像に対して位置合わせされた規則的グリッド上でワーピングパラメータ（並進、回転等）を計算する。
2. 前記規則的グリッド上での補間によって各像点のワーピングパラメータを求める。湾曲歪みを補正する処理が合成画像に対して及ぼしうる効果が図11に図示されている。図11の観察から、湾曲歪みが補正されたSEM画像から得られた合成画像（図11B）の品質は、湾曲歪みが補正されない場合に得られる合成画像（図11A）の品質よりも優れていることは明らかである。
【００５６】
上述した湾曲歪みの補正方法は例示であり、他の方法を用いることも可能であることに留意して欲しい。
【００５７】
C. 互いにレジストレーションされた複数の3D画像の処理
任意に配向（たとえば回転、伸縮、シフト）する3D画像を3D画像に対して正確にレジストレーションすることに依拠する他の応用についてここで説明する。結果として得られた合成情報は、多くの複雑な実験が有孔性材料について行われた後、ある範囲の条件下で同一の有孔性材料についての物理的特性の測定を可能にすることによって、従来のマイクロCTイメージングの能力を大幅に拡張する。
【００５８】
かなりの長期間にわたって得られた一連の3D画像が3Dでレジストレーション可能である場合、流体相分布は、顕著な劣化の後又は長い平衡時間の後に、3Dにおいて孔のスケールで直接的にマッピングされて良い。よって全範囲の所望の物理特性が画像データに基づいてシミュレーションされて良い。これらには、相対透水率（有孔性材料内部の複数の流体相の透水率）、複数の流体相が集まった岩石の抵抗力、弾性応答特性、NMR応答、及び動的毛管圧力（たとえば非特許文献3を参照のこと）が含まれる。流体分布は、孔ネットワークモデルの結果と直接比較することが可能で、かつ孔スケールのモデル化手法の精度を上げるための校正ツールとして用いられて良い。
【００５９】
開示された方法は概して、孔構造、鉱物相構造、及び3次元有孔性媒質の孔空間内部での複数の流体相の分布の変化の視覚化並びに定量化に関する。有孔性媒体の変化の直接的な定量化から、様々な飽和状態及び様々な表面化学状態下での有孔性材料のある範囲の物理特性が推定されて良い。これらの特性には、様々な流体相の相対透水率、抵抗力、誘電応答、NMR応答、及び弾性率が含まれる。材料には、油層のコア材料、土壌、粒状材料、及び他の材料が含まれる。
【００６０】
開示された方法は、現実的な条件下での多相流、岩石の機構、及び油層岩石物理学上の特性の理解の必須要件を解決する。たとえば、流体の吸水特性はコアの未精製油の経年劣化後に変化することは既知である（たとえば非特許文献7,8を参照のこと）。他の測定もまた、油層温度及び岩圧での実験を行えることが求められる。提案された方法は、説明したアルゴリズムによって実装されたソフトウエアによる、試料の撮像、実験室中での前記試料のある範囲の複雑な実験の実行とそれに続く前記試料の再撮像、及び、前記実験前後で得られた3D画像のレジストレーションを可能にする。
【００６１】
レジストレーションされた複数の3D画像に係るデータの合成方法
この方法の手順が図12に図示されている。当該方法は、関心対象である有孔性材料の3D画像を取得する手順1202で開始される。この撮像は最初の参照状態Aを定義する。一般的には、この状態で研究される試料は、乾燥状態（図13）又は飽和流体（図14）である。画像データは手順1204においてメモリ内に記憶される。
【００６２】
続く手順1206は、孔構造、鉱物相構造、若しくは孔への流体分布を構造的又は科学的に変化させることが可能な有孔性材料上での1つ以上の実験の実行を含む。実験例には以下が含まれて良い。
− 有孔性材料についての流体変位の研究。該研究は以下のようなものであって良い。
・ 非混和性流体の変位
・ 流体の侵襲
・ 界面活性剤、微生物、ポリマー、ゲル、又はコロイドによるフラッディング
・ 混和性流体の変位
− 反応による効果。該効果は以下のようなものであって良い
・ 反応流
・ 溶媒によるフラッディング
・ 堆積及び反応による汚損
・ 生物付着
− 力学的効果。該効果には以下のようなものがある。
・ グレインの鉱物学的特性/破壊
・ 粘土の膨張、付着、及び分離
・ 流体の音響学的シミュレーション
− 孔及び鉱物相の構造への流体の効果
・ 流体の交換
・ 微粒子群、ポリマー、及び/又は水性溶媒によって孔を塞ぐこと
・ 粘土のマイグレーション（天然）
・ コロイドの付着/分離及び汚損
・ 粘土の膨張
− 濡れ性分布。該濡れ性分布には以下が含まれる。
・ 流体分布を孔のサイズ及び構造に相関させること
・ 濡れ性及び孔への流体の分布に対する未精製状態での経年劣化の影響
・ 液体/液体及び液体/固体界面の表面化学特性の測定（たとえば接触角、流体膜の存在）
当該方法は、手順1208において、1つ以上の実験を行った後に前記試料材料を再撮像する手順、及び、手順1204において、前記の再撮像された画像を記憶する手順をさらに有する。前記試料は複数回再撮像されて良い。前記試料による実験を行う代わりに、手順1206は、より高い解像度での試料のさらに小さく分割された部分の取得に係るものでも良い。その結果、実験前後で取得された画像又は様々な撮像解像度で得られた画像のいずれかを表す1対の画像が生成される。よって2組のデータは、3D画像のボクセルレジストレーションを実行することによって得られた画像のうちの2つから得られた3Dデータを位置合わせして、かつ重ねることによって、合成画像情報を得るのに用いられて良い。そのような位置合わせ及び合成は、続いて別な画像の対の間で行われても良い。2つ以上の画像の情報を合成する処理は、両画像からの情報に基づいて、前記試料の特定領域の物理特性又は化学特性を決定する手順を有して良い。当該処理はまた、孔の構造、鉱物相構造、又は孔への流体分布の変化を定量化する手順をも有して良い。
【００６３】
ここで図13と図14を参照すると、図は、流体取り込み後の有孔性岩における空気/水相を図示している。3Dレジストレーションが実行されることで、流体の孔スケール分布を3Dで直接的に視覚化及び定量化することが可能となる。特に図14は、乾燥した岩へ濡れ性流体が取り込まれる実験を行った後、図13から得られた岩石の流体分布を図示している。
【００６４】
図15は、一対の3D画像データファイルの画像レジストレーションから得られた流体分布を表す3Dビデオからのスナップショットを図示している。
【００６５】
試料材料の再撮像は、各実験が行われた後又は一連の実験が行われた後、何度も繰り返されて良い。繰り返しになるが画像データは、コンピュータのメモリ又はデータ処理を取り扱う他の電子システム（通常はハードドライブ）に記憶される。よって、異なる画像を取得する間の期間に制限はない。よって2つ以上の記憶された画像のボクセル−ボクセルレジストレーションは、任意の数の画像から得られる情報を結合させることで、全実験値に追加することを可能にする。
【００６６】
任意の期間にわたって行われる任意の数の実験から合成画像を生成するために提案された方法は、従来技術を有する2つの実験方法とは異なる。非特許文献9に代表される前者の方法では、全ての実験はビームライン中で行われる。このことはつまり、実験の性質及び期間は、ビームラインの利用可能性、空間についての考慮、並びに、温度及び機械的安定性の要求によって制限される。後者の方法（たとえば非特許文献10を参照のこと）では、試料が撮像され、続いて実験の修正を行うために移され、その後再撮像のためにビームラインに戻される。しかし異なる画像のレジストレーションを行うための努力は行われない。この場合では、合成画像の情報は直接得られず、唯一の結果は空間的な平均である。
【００６７】
最終手順1212は、手順1210においてレジストレーションされた3D画像から得られた合成情報を用いることで、孔の構造、鉱物相構造、又は孔への流体の分布の変化を定量化する手順を有する。これらの変化の有孔性材料の物理特性に対する示唆が調査されても良い。たとえば濡れ性分布の孔スケール分布の直接的な定量化から、様々な飽和状態及び様々な表面化学条件下でのこれらの有孔性材料のある範囲の物理特性の直接的なシミュレーションが検査されて良い。関心特性には典型的には、様々な流体相の相対透水性、抵抗力、誘電応答、及び弾性が含まれる。
【００６８】
図16A及び図16Bは、孔のサイズ分布と空気/水流体の含有量との間の相関を示している。3Dレジストレーションにより、孔スケールでの流体の分布を3Dにおいて直接的に可視化及び定量化することが可能となる（図16A）。図16Bは、濡れ性を有しない（空気）相が3Dにおいて存在する場所の分布を示している。この場合では、系は本来水で濡れているので、残りの空気相は基本的に大きな穴にとどまり、かつ水相は小さな孔に集中する。
【００６９】
D. 典型的な応用
孔の構造、鉱物相の構造、若しくは孔への流体の分布の起こりうる構造変化又は化学変化を定量化することが可能であることが重要となる、ある範囲の学問分野にわたる有孔性材料に関連する多数の実験がここで論じられる。
【００７０】
1. 孔のスケールでの有孔性材料の流体変位の記述。有孔性材料の撮像の第1用途の中には、様々な濡れ性条件及び飽和状態の下で記載され、かつ再撮像されるものがある。複数の相が有孔性材料内部で撮像されて良い。流体変位特性に関連する実験には以下が含まれて良い。
【００７１】
a. コア解析における非混和性流体変位（上流の油とガス/地下水及び水門学）
・ 複数の飽和状態での排水及び吸水フラッディング
・ 様々な濡れ性状態（油が濡れる、水が濡れる、混合物の濡れが大きい、混合物の濡れが大小さい）下での排水及び吸水
・ たとえばフラッディング中及びフラッディング後における、未精製油の経年劣化が、炭化水素、ブライン、ガスなどに及ぼす効果の探索
・ ある範囲の流体分布の探索。
− 油/炭化水素
− 水/ブライン
− CO2
− 泥
− 入り込む流体
− ポリマー
− コロイド
− 分散剤
− 孔へ入り込む流体（たとえばポリマー、ゲル、微生物種）
・ 意図した充填の試験
− 孔構造における抽出を最大にするための最適なブライン
− 様々な油相/ブラインの混合物での孔スケール分布の試験
− 改良/改善された油の回復戦略（たとえばWAG注入及びコロイドシミュレーション）の効率
・ コア材料への泥の侵襲（損傷の生成）の撮像
b. 混和性流体変位
・ 混和性フラッディング(EOR)
・ 地下水の浄化
− 孔スケールでの毒物又は放射性物質の広がりの分布の定量化
− 汚染物（化学的、水門学的）の試験方法
2. 孔及び鉱物相の構造への流体の効果：孔の構造、鉱物の構造、及び有孔性材料の濡れ性に対する流体の交換の効果（たとえば孔の幾何学形状に対するブラインの希釈効果）
a) 粘土のマイグレーション（天然）
b) コロイドの付着/分離及び汚損
c) 粘土の膨張
3. 濡れ性の分布。多数の流体相によって占められる有孔性材料の流体/流体界面及び流体/固体界面の濡れ性状態が評価される。これらには以下が含まれる。
a. 流体分布を孔のサイズ及び構造に相関させる
b. 液体/液体界面及び液体/固体界面の表面化学特性の測定（たとえば接触角、流体膜の存在）
4. 反応の効果。有孔性材料の孔及び鉱物相の構造への反応種の効果を直接孔のスケールで調査する。流体の濡れ性についての研究は複数回の手順で行われても良い。
a. 反応流
− 超臨界CO2
− 酸のフラッディング
− 溶媒に基づくフラッディング
b. 堆積及び反応による汚損
c. 生物汚損
5. 力学的効果。有孔性材料の孔及び鉱物相の構造への力学的な応力と歪みの効果。これらには以下が含まれる。
− グレイン鉱物学的効果/破壊
− 粘土の膨張、付着、及び分離
− 流体の音響試験。炭化水素の回復を改善するための音響波の試験的利用。
【００７２】
E. 画像レジストレーションのワークフロー
以降で説明するレジストレーションのワークフローは、2D-3D画像レジストレーションにも3D-3D画像レジストレーションにも適用可能であることに留意して欲しい。定義を多少変更させれば同一のレジストレーションのワークフローは2D-2D画像レジストレーションにも適用可能である。
【００７３】
定義
Tを1つの独立の画像全体を表すものとし、Ti,j,kを指数(i,j,k)での前記独立の画像の強度を表すものとする。Fを他の独立の画像全体を表すものとし、同様にFi,j,kを指数(i,j,k)での前記他の独立の画像の強度を表すものとする。さらに、I(n)を、dn倍のダウンサンプリングを行うことによって独立の画像Iから生成される独立の画像であるとする。T(0)⇔T及びF(0)⇔Fとなるように常にd0=0と定義する。独立の画像Tは動く独立の画像で、かつFは固定された独立の画像である。Tが動くことで、Tは固定された独立の画像に対してレジストレーションされる。たとえばTは、3DマイクロCT画像であるFに対してレジストレーションされる2DのSEM画像であって良い。Tが2D画像である場合、Tは3次元空間内の一面であると単純に解される。つまり第3の指数は、k=0に固定される。あるいは別な例として、Tは、3Dの乾燥状態のマイクロCT画像であるFに対してレジストレーションされる3Dの湿った状態のマイクロCT画像であって良い。独立の画像FとTはまた一対の2D画像であっても良い（たとえば図17に図示されたSEMの薄い部分の画像の一部分174）。またわずかな回転によって変化させることで、後述の方法は2D-2Dレジストレーションにも適用される。
【００７４】
Ωt⊂R3であるとき、t(n):Ωt→Rを独立の画像T(n)の補間された画像T(n)であるとする。同様に、Ωf⊂R3であるとき、f(n):Ωf→Rを独立の画像T(n)の補間された画像F(n)であるとする。
【００７５】
Ψ(n):Rm(n)×R3→R3を空間変換と定義する。Ψ(n)は、一組の変換パラメータφ∈Rm(n)が与えられると、3Dの点を3Dの点へ写像変換する。
【００７６】
M(n):Rm(n)→Rを、2つの補間された画像t(n)とf(n)との間でのレジストレーションの品質を示す距離と定義する。またM(n)の値が小さければ小さいほどレジストレーションは良好で、M(n)の値が大きければ大きいほどレジストレーションは不十分である。
【００７７】
初期の独立の画像t(0)とf(0)が与えられると、2つの画像をレジストレーションさせる変換パラメータφを求める必要がある。これは最適化問題に再定式化することができる。
【００７８】
レジストレーション最適化問題：空間変換Ψ(0)が与えられると、距離M(0)及び画像t(0)とf(0)によって、φ*∈Rm(0)が求められる。このとき、
【数１】
が満足される。
【００７９】
一般に、M(0)は多くの極小値を有する複雑な関数で、かつ（特にフル解像度の画像について）評価するのには計算上手間がかかる。従って、局所的な数値最適化法（シンプレックス法、パウエル法、勾配降下法等）単独ではφ*（の近似値）を巧く求めることができない。以降では、レジストレーション最適化問題を解決することが可能な多解像度多出発点の大域最適化法が与えられる。
【００８０】
一般的方法
1. ダウンサンプリング因子d0,…dN、空間変換Ψ(0),…Ψ(N)、距離M(0),…M(N)、及び最善の変換パラメータの数J(0),…J(N)を選ぶ。たとえば以降のこれらの変数の選択を参照のこと。
2. φ(N)⊂Rm(N)を、Rm(N)から得られた有限個のサンプリング点と定義し、かつP(N)を、組φ(N)中の要素数とする。各φP(N)∈Φ(N)、p=1,…P(N)について、mp(N)=M(N)(φp(N))を評価する。以降のたとえばΦ(N)を参照のこと。
3. 組Φ(N)（ハット）={φj:φj=φpj(N)、j=1,…J(N)}を生成する。ここでpjは、mpj(N)がすべてのmp(N)、p=1,…P(N)のうちでj番目に小さな値であるように定義される。
4. n=Nと設定する。
5. n=n-1と設定する。n≧0の場合には6へ続き、そうでない場合には9へ進む。
6. j=1,…J(n+1)の各々について、開始点φj(n+1)∈Φ(n+1)（ハット）を有する局所的数値最適化アルゴリズム（たとえばパウエル法）を用いてM(n)の極小値を計算する。極小値を計算するための変換パラメータをφj(n)（ハット）とし、かつこれらの極小値での距離をmj(n)（ハット）=M(n)(φj(n)（ハット）)とする。
7. 変換パラメータの組Φ(n)（ハット）={φj:φj=φpj(n)（ハット）、j=1,…J(n)}を生成する。ここで、pjは、mpj(N)（ハット）がすべてのmp(N)（ハット）、p=1,…P(N)のうちでj番目に小さな値であるように定義される。
8. 5へ進む。
9. 次式を、t(0)の画像をf(0)の画像へレジストレーションするのに最適なパラメータの組と設定する。
【００８１】
【数２】
10. 前記最適なφ（ハット）の変換パラメータを用いることによって、レジストレーションされた状態である独立の画像の新たな対を生成する。
【００８２】
手順1は、ダウンサンプリング因子、許容される空間変換、及び所与の組の変換パラメータで実現されるレジストレーションの品質評価に用いられる距離の種類についての詳細を定義する。典型的には、ダウンサンプリング因子は2の乗数となるように選ばれる。なおd0=1で、n0<n1であればdn0≦dn1である。例の結果、d0=1,d1=2, d2=4,d3=8,
d4=16,及びd5=16である。最大のダウンサンプリング因子は、特徴部が最低解像度の画像となるのが防止されるように選ばれなければならない。
【００８３】
空間変換Ψ(0),Ψ(1),…Ψ(N)は、どのようにしてt(n)の領域Ωtからの座標が、f(n)の領域Ωfへレジストレーション/重ね合わせされるのかを定義する。一の可能性は、すべての空間変換が等しくなるように、Ψ(0)⇔Ψ(1)⇔・・・⇔Ψ(N)⇔Ψ（オーバーバー）と選ぶことである。これは、図中の例をレジストレーションする場合には正しく、かつΨ（オーバーバー）は3D相似変換となるように選ばれた（Ψ（オーバーバー）:R7→R、7の自由度は、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータを有する）。
【００８４】
選択可能な実際の距離の指標には様々なものが存在し、繰り返しになるが、妥当な選択とは、各探索レベルn,M(n)（オーバーバー）の距離について同一の数式を有することである。つまり探索レベル間での唯一の違いは、距離が対応するダウンサンプリングされた画像対を用いて評価されることである。有用なコストの指標には、相関係数、相関比、及び（規格化された）相互の情報が含まれる。図中のレジストレーションされた例については、M(n)（オーバーバー）が負の相関係数として選ばれた、
【数３】
ここでt(n)（オーバーバー）は平均で、σt(n)は試料座標xlでのt(n)画像の強度の標準偏差で、f(n)は平均で、かつσf(n)は試料座標Ψ(φ,xl)でのf(n)画像の強度の標準偏差であり、試料の空間座標は画像領域の交差部分xl∈(Ωt∩(Ψ-1(φ)（オーバーバー）・Ωf))、から取得される。ここでΨ-1（オーバーバー）はΨ（オーバーバー）の逆関数で、Ψ（オーバーバー）（Ψ-1（オーバーバー）(x)）=x,∀x∈R3となる。
【００８５】
J(n)の値は、探索レベルnから探索レベルn-1へ渡される最善のランク付けがなされた変換パラメータの数である。典型的には、J(0)=1で、かつn0<n1のときにはJ(n0)≦J(n1)となる。J(n)の結果の一例は、J(0)=1, J(1)=1, J(2)=1,
J(3)=4, J(4)=32, J(5)=64である。
【００８６】
手順2では、目標は、所望の変換パラメータφ*の組に「近い」変換パラメータの組をすくなくとも1組求めることである。この手順では、包括的な型の探索が、より高い解像度での所望のレジストレーションを実現すると思われる変換パラメータの組を発見するために実行されて良い。典型的には、Φ(N)の組は、変換パラメータΦ(N)={φα:Zm(n),α=0,…A}の正方形グリッドとして定義される。ここでα=[α0,α1,…αm(N)]Tは多重添字である。φαは次式に示す行列で表される。
【００８７】
【数４】
ここでφmin∈Rm(N)は変換パラメータ成分の下限で、δ∈Rm(N)は変換パラメータ空間内のグリッドのステップサイズである。最小の並進範囲は、Ωf∩Ψ(N)(φ)・Ωtが全てのφ∈Φ(N)に対して空集合とならないことを保証するように選ばれて良い。ここで典型的な並進ステップサイズはダウンサンプリングされたボクセルサイズである。回転値は-180°から180°の範囲で、ステップサイズは約3°である。手順2は基本的には、距離の大域最小の「捕獲半径」の範囲内にある1組の変換パラメータを求めるために低解像度の画像に適用される大域最適化法である。この段落は、「包括的な探索」による大域最適化法について説明している。しかし最適レジストレーションパラメータに「近い」変換パラメータ（又はその複数の組）を発見するため、他の大域法も手順2において適用されて良い。
【００８８】
手順3は、最高のレジストレーションから最低のレジストレーションまでの変換パラメータの組を順序づける。繰り返しにより、手順4から手順8は、順次より高解像度の画像へ局所最適化を実行することによって、手順2で求められた変換パラメータの組を徐々に改善する。
【００８９】
手順9は最適なレジストレーションを与える変換パラメータを選ぶ。手順10は、画像の再サンプリングを（その画像の最適変換パラメータを用いて）行うことで、後続の解析を行うために直接比較することが可能な一対のレジストレーションされた独立の画像を生成する。
【００９０】
アルゴリズム
上述のレジストレーション法を実装するためのアルゴリズムが図17に図示されている。図17では、3D画像1702と2D画像1708との間でのレジストレーションが実現される。最初に得られた画像1702は、フィルタリング、マスキング等によって、プロセッサ3605Aによって前処理されることで、レジストレーションの準備が整った画像1704が得られる。同様に、個々の画像1706もまた前処理されることで、2Dの操作された画像が得られる。特に、2Dの再分割された部分1706の前処理は、湾曲歪みの補正及び分割画像を一つにまとめることで1つの画像1710を得る手順を有して良い。
【００９１】
図17のアルゴリズムの手順1では、2D画像1710を3D画像1704に対して位置合わせするため、相似変換がプロセッサ3605Aと3605Bによって選ばれる。7の自由度にはこの変換も含まれる。具体的には、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータである。負の値を有する相関係数は、位置合わせの品質を評価するのに用いられる。-1の値は完全な相関を意味する。
【００９２】
手順2では、一対の低解像度の画像が、プロセッサ3605Aと3605Bによって用いられることで、画像1704と画像1705との間での位置合わせについての一組の初期推定が得られる。手順2の実装は、タスクの並列化によってその手順の実行時間を減少させるという利点を得ることができる。この用途におけるタスクの並列化の実装についてのさらなる詳細は本明細書においてさらに与えられる。手順3では、局所最適化用の開始点として用いられる一群の最善の変換パラメータが選ばれる。当該処理のこの段階では、画像1704と画像1710とは最低解像度の状態で大まかな位置合わせがなされている。
【００９３】
手順4では、多数の探索に関してのループカウンタが、プロセッサ3605Aと3605Bによって初期化される。続く手順5では、ループカウンタは、プロセッサ3605Aと3605Bによってすぐに1つずつ減ってゆく。ループカウンタが0以上の値を保持する場合、当該方法は手順6へ進む。手順6では、プロセッサ3605Aと3605Bは、現時点での高解像度画像で局所最適化を用いることによって、過去の探索レベルn+1からの位置合わせを精緻化する。手順6の実装は、低解像度の画像についてのタスクの並列化を行うという利点を得ることができる。前記タスクの並列化では、3605Aと3605Bがそれぞれ異なる複数の開始点の部分集合を処理する。より高解像度の画像の対については、画像に必要なRAMが1つのコンピュータプロセッサの物理RAMを超える場合、データ並列化を実装することができる。前記データ並列化では、3605Aと3605Bはそれぞれ、各プロセッサからアクセス可能なメモリ内に記憶された各プロセッサに対応する部分像に関する指標の各独立した一部分を計算する。係る実装についてのさらなる詳細は本明細書においてさらに供される。手順7では、プロセッサ3605Aは、次の探索レベルの開始点として用いられる現在の探索レベルから最善の変換パラメータの群を選ぶ。手順8では、そのループは手順5へ戻ることによって続けられる。
【００９４】
nが負の場合には他のルートがとられる。この場合では、手順5の後に、プロセッサ3605Aと3605Bは手順9へ進む。そこでプロセッサ3605Aと3605Bは探索レベルn=0から最善の変換パラメータを選ぶ。これらは、フル解像度での最低指標値を供するパラメータである。これらのパラメータを用いることで、手順10では、プロセッサ3605Aと3605Bは、3D断層画像を再サンプリングすることで、2D画像に対応する薄片を生成する。この結果、同一サイズである2つの位置合わせされた2D画像が生成される。よってこれらの画像は、後続の定性的解析及び/又は定量的解析を行うのに直接用いられて良い。
【００９５】
図17の手順1から手順10は、図1で開示された画像処理方法の高次表現に係る手順108と指称されるレジストレーション方法を実質的に規定する。同様に図18は、図12で開示された画像処理方法の高次のダイアグラムに係る手順1210を示している。
【００９６】
図18に示されたレジストレーション方法は、2つの3D画像1802と1808にも適用可能である。画像1802は乾燥状態のコア試料の3DマイクロCT画像であって良い。対照的に、画像1808は、特別の実験を行った後の同一試料の3DマイクロCT画像であって良い。前記特別の実験にはたとえば、特別の流体−たとえば水、石油、又はガス−へコアを導入する手順が含まれて良い。画像1802と1808のいずれも、たとえば単純な閾値設定を利用したフィルタリング及びマスキングによって前処理される。前記前処理の結果、画像1804と1810がそれぞれ得られる。繰り返しになるが、手順1では、プロセッサ3605Aは、2つの断層画像1804と1810の位置合わせを行うため、相似変換を用いる。繰り返しになるが7の自由度がプロセッサ3605Aによって選ばれる。前記7の自由度は、3の並進パラメータ、3の回転パラメータ、及び1の伸縮パラメータを有する。負の重なり割合が用いられる。この負の重なり割合は、両画像での同一のマスキングされていない強度を有する座標の割合を単純に計数するものである。-100%の値は両画像においてマスキングされていない領域が完全に重なっていることを意味する。
【００９７】
手順2から手順9は概して図17に図示された手順と同一である。繰り返しになるが手順10は、最終手順であり、かつ、乾燥状態の3D断層画像1804に対応する領域を生成するため、湿った状態の3D断層画像1810を再サンプリングするプロセッサ3605Aと3605Bを有する。その結果、一対の3D画像1814と1816が、プロセッサ3605Aと3605Bによって得られる。一対の3D画像1814と1816は、同一サイズであって、かつ定量的な解析目的のために合成されてよい。
【００９８】
図19は、図17及び図18の手順2の様々な下位手順をより詳細に図示している。手順2は、7次元の変換パラメータ空間において規則的なグリッドをプロセッサ3605Aと3605B によって生成する手順を有する手順2aで開始される。典型的には並進パラメータのステップサイズは最低解像度の画像のボクセルサイズ長である。z軸の回転パラメータのステップサイズは2°である。x軸とy軸の回転パラメータ及び伸縮パラメータは（それぞれ0°、0°、及び1で）固定されて良い。一旦手順2aから規則的なグリッドが生成されると、プロセッサ3605Aは、kを7Dパラメータグリッドにおける点の数とすることによって、パラメータの初期化を行う手順2bへ進む。続く手順2cでは、パラメータkは1だけ減少し、かつプロセッサ3605Aと3605Bによって、新たなパラメータ値が0と比較される。
【００９９】
減少したkが0以上である場合、プロセッサ3605Aと3605Bは、手順2cから手順2dへ進む。手順2dでは、プロセッサ3605Aと3605Bはそれぞれ、7Dグリッド点の数kの変換パラメータの指標を評価する。前記指標の値がこれまでに得られた最善の指標値よりも低い場合、プロセッサ3605Aと3605Bは、最低の指標値を与える変換パラメータを変換パラメータ数kに置換する。手順2dの結果は、プロセッサ3605Aと3605Bを手順2cへ戻すためである。あるいはその代わりに、パラメータkの値が、手順2cで減少したように、負である場合、プロセッサ3605Aと3605Bは手順2cから手順2eへ進む。手順2eは手順2の全体を終了させる。パラメータkが負の値をとるということは、7Dグリッドの包括的な探索が完了し、かつ指標となる大域最小値の捕獲半径の範囲内に属する少なくとも一組の変換パラメータが特定されたことを意味する。ここでタスクの並列化は、各異なるプロセッサ−たとえば3605A、3605B等−間のkのグリッド点を分割することによって実装されて良い。各プロセッサは、グリッド点の部分集合についての指標を評価する。最善のパラメータの群は、各プロセッサから収集され、かつ照合されることで、最善の変換パラメータからなる最終群が生成される。
【０１００】
並列実装
上述の大域最小化法を計算機上で実行可能にするため、高性能非均質メモリアクセス(NUMA)プラットフォーム上のソフトウエアで上述の手順を実行されるときに使用される2つの並列化法が存在する。第1の方法はタスクの並列化法である。前記タスクの並列化法では、個々のタスクは計算ユニット上で独立に実行される。第2の方法はデータの並列化法である。前記データの並列化法では、各計算ユニットは独立の画像対データの部分集合だけしか含んでいない。
【０１０１】
タスクの並列化法は、上述のワークフローの手順2において距離の値mk(N)を計算するときに大いなる利点となる。低温では、各計算ユニットは、他の計算ユニットとは独立してmk(N)の値の部分集合を評価することができる。全ての計算ユニットが、それらのmk(N)の値の部分集合を評価したとき、結果として収集段階では、手順3で実行されるように、最善の変換パラメータを順位付けするため、各計算ユニットの各々からmk(N)の値の部分集合が集められる。同様に、各計算ユニットが極小値Φj(n)（チルダ）を計算するため、他の計算ユニットとは独立にして繰り返し最適化の一部を実行することができる手順6で利点を得るためにタスクの並列化が用いられて良い。繰り返しになるが、結果として実現される収集段階は、手順7で実行されるように、各計算ユニットからのmk(N)の値の部分集合を収集して順位付けするのに用いられる。タスク並列化は低解像度の画像（たとえばt(N), t(N-1), f(N),及びf(N-1)）にとって有利である。なぜなら基礎となる独立の画像対の全データは計算ユニットの局所RAM内に記憶することができるからである。
【０１０２】
高解像度の画像（たとえばt(0), t(1), f(0),及びf(1)）にとっては、データの並列化法が好ましい。この手法では、各計算ユニットは独立の画像のデータの部分集合しか有していない。よって手順7は、マスターワーカー並列化によって実行される。マスター計算ユニットは、繰り返し最適化アルゴリズムを制御する。前記繰り返し最適化は、極小値を決定するため、指標となる距離の多重評価を必要とする。指標となる距離の評価は、分配して、かつ並列に行われる。各ワーカー計算ユニットはまた、該ワーカー計算ユニットの局所RAM内に保持された画像データの部分集合に基づいて、全体の指標となる距離の一部を計算することができる。
【０１０３】
E. さらなる応用例
X線マイクロ断層撮像(CT)は、ミクロンスケールでの孔の構造の詳細な3D画像を生成する能力を有する。X線断層撮像は、検出されたX線の減衰を介して位相情報を供する。このデータについては詳述しないが、様々な2D顕微鏡手法によって供されるものとして有用である。たとえば光学顕微鏡及び走査型電子顕微鏡は、試料の詳細な鉱物学的内容の情報を供することができる。走査型音響顕微鏡は、鉱物の音響（弾性）特性のマッピングを可能にする。2D表面又は薄い部分での実行が可能な他の顕微鏡手法は、研究対象であるコア材料（試料）の材料特性の理解にとって重要な他の情報（たとえば鉱物学的知見、表面特性、局所的表面相互作用等）を与えることができる。一の欠点は、得られたデータが2Dであることである。2D方法によって利用可能な情報が3Dにおいて改善された情報として得られることは、大きな利点である。ここで述べている方法は、2D画像から得られた情報を3D断層画像データの組に定量的に伝播することを可能にする。
【０１０４】
上述のレジストレーション方法はこの伝播処理の一部である。図20は2D画像2010を3D画像2012へレジストレーションした結果を図示している。特に、当該方法において定義されている全ての並進、回転、ワープ、及び伸縮変換が行われた後、2D顕微鏡で得られた薄い部分2010は試料2012内部でレジストレーションされることで、CT画像2012から得られた3Dデータ組への、2D画像2010から得られたSEMデータの結合が可能となる。
【０１０５】
マイクロ有孔性岩石の試料の3D断層画像への高解像度2DのSEM画像から得られた画像データのレジストレーションの例が図21に図示されている。2DのSEM画像2112から得られたマイクロ有孔性データを3D断層画像2110へ結合することによって得られた合成データは、2DのSEM試料での解像度、及び3D断層撮像によって可能な解像度よりも高い解像度で、材料の部位及び特性を識別することを可能にする。一般的には、断層撮像手法は、本質的な解像度の限界（たとえば1μm）を有する。顕微鏡はnmスケールでの探索を可能にする。よって、画像レジストレーションによって2つの画像を合成することで、最大1000倍の解像度のオーダーが利用可能となる。それに加えて、2つの手法の解像度の大きな差異は、3D画像の真の解像度を測定するのに用いることができる。
【０１０６】
図22は、岩石試料の3D断層画像から得られた2D薄片2210と、SEMから得られた2D薄片との直接比較を図示している。断層画像2210は2.1μmのボクセルサイズのX線CT撮像によって得られているのに対し、SEM薄片2212は0.25μmの解像度で得られている。2つの画像2210と2212との定量的な一致が図23に図示されている。具体的には、図23は、高解像度の2D画像2212から得られた（特定の領域の一部に係る）有孔率データを3D薄片2210から得られた減衰データへマッピングするグラフを図示している。グラフ2314から分かるように、減衰データと有効率データとの間には良好な相関が得られた。このため、高解像度2D顕微鏡像から得られた情報（この場合は有孔率）を、3D全体へ伝播させることが可能となる。これは、3D画像の減衰データに基づくものであり、かつグラフ2314によって図示された直接的かつ定量的にマッピングされた相関を用いている。
【０１０７】
2DのSEM画像から高解像度で得られる他の特性もまた、レジストレーションされた3D画像と定量的に合成されて良い。図24A及び24Bに図示された例では、SEM画像2410は0.25μmのピクセル解像度で得られる。画像2414は、2.5μmのボクセルサイズで得られたX線CT画像からのレジストレーションされた薄片である。この画像でのX線減衰は、画像2410内の孔のサイズデータに相関付けられて良い。画像2412は、画像2410から得られた画像データ上での局所的な孔のサイズ分布のマッピングを表す。画像2412中でのグレースケールは、2410から得られた孔のサイズに対応する。明るい領域は大きな項のサイズに対応する。
【０１０８】
図24Bは、局所的な孔のサイズ分布（y軸）とX線減衰からのマイクロ有孔率マッピング（x軸）との相関のプロットを図示している。この相関により、3D画像全体への孔のサイズをマッピングすることが可能となる。
【０１０９】
図25は、2D顕微鏡画像から3D断層画像への定量的な画像レジストレーションの別例を図示している。この例では、元素マッピング（画像2512）から得られた鉱物相の情報が、画像2510によって表される3D画像から得られたレジストレーションされた薄片で得られた減衰データに関連づけられている。2D画像2512に示された詳細な鉱物相の分布は、岩石試料内でのドロマイト及び硬石膏の存在を示している。これらの相は、3D画像2510では明瞭に差異化するのが難しい。両画像間でのレジストレーションを利用することによって、定量的マッピングが、3D画像全体での鉱物相の情報を識別する能力を改善させた。
【０１１０】
ここで、3D-3Dレジストレーションに係るさらなる応用について説明する。係るレジストレーションによって、1つ以上の実験が試料上で行われた後に、その試料中に生じた変化を定量的に追跡することが可能となる。3D-3Dレジストレーションを含む典型的な応用は通常、試料−特別な状態Aにある−の撮像で開始される。この最初の画像が得られた後、試料は撮像装置から取り外され、かつ、様々な実験−たとえば圧密、圧縮、分解、流体変位−が様々な条件下で行われる。係る実験によって試料の状態（たとえばマイクロ構造、鉱物学的状態、流体飽和）が変化する。続いて、この時点で新たな状態Bにある試料は、撮像システムへ再導入されて再撮像される。続いて新たに得られた画像が先に得られた画像に対してレジストレーションされる。2つの画像から得られたデータは合成され、生じた変化は定量化することができる。
【０１１１】
以下の例は、一組の実験中での分解と粒子運動の視覚化及び定量化を有する方法の応用、並びに様々な化学条件下での複数の流体飽和状態の撮像を表している。また複数のスケールでの3D画像の位置合わせ例も与えられている。
【０１１２】
図26及び図27は、一連の実験中に試料体積内部の物理構造において生じた変化を定量的に追跡する実験を表している。係る塩化は鉱物に固有であって良い。岩石試料2610は、実験室で受け取られたような本来の状態である。薄片2616は、試料2610が本来の状態であるときに、その試料2610の3D画像内部から得られたものである。続いて試料2610には、中間状態2612によって示されているように、1つ以上の実験が行われる。図26の場合、試料2610は特定の組成を有する流体に浸される。その流体は、特定の期間中高温かつ高圧に維持される。特別なセル内で維持されたそのような条件は、試料材料の分解及び沈降を起こす。その結果、長期間にわたって試料のミクロ構造が変化する。そのような変化を評価するため、試料は新たな状態2614で再撮像される。3D画像内部からの薄片2618もまた図26に図示されている。2616及び2618のより詳細な拡大図が、2720と2722にそれぞれ図示されている。薄片2720と2722とを比較すると、試料の形態において視覚上明確な差異が存在する。前記差異は特定の実験条件によって引き起こされたものである。薄片2720と2722に係る2つの大きな3D画像を合成する結果、多数のデータが、遷移期間中に試料体積全体にわたって生じたミクロな変化を捉えることができた。例の小さな部分集合は、2616から2618への遷移中に発生した、孔構造、グレイン構造、又は、流体透水率の変化の列挙を有して良い。これらの局所スケールの差異/変化は、当該方法に含まれる直接レジストレーションによってしか識別できない。
【０１１３】
図28は、カルボン酸への曝露後における岩石の分解/沈降特性の変化を研究するための様々な実験を図示している。薄片2810は初期状態における岩石試料の本来の3D画像を表している。薄片2812は、長期間にわたってカルボン酸へ曝露した後の同一岩石試料からとったレジストレーションされた3D画像から得られた薄片を表す。前記曝露の前後で2つの3D画像のレジストレーションを行うことによって、岩石試料中で生じた形態変化を識別することができる。この特別な例における変化は顕著でかつ明らかである。たとえば、本来の画像2810において視認可能な部位2814は、酸への曝露後に撮られた画像から得られた、同一で定量的にレジストレーションされた薄片2812から消滅した。よって、2つの画像間での差異を示す合成された3Dデータの組は、材料構造の変化を定量化することを可能にする（たとえば孔の構造、孔の形状、グレインの形状、ミクロ有孔率）。物理的特性の変化もまた定量化することができる。たとえば分解及び再沈降による試料の流動特性の変化は、合成された3Dデータ組を用いることによって探索可能である。地質材料の分解を調査するための様々な方法について検討するため、さらなる実験が様々な条件下で行われても良い。
【０１１４】
図29は、オイル、ガス、及び/又は水をコア試料全体へフラッディング（循環）させる間に生じる様々な状態によって生成される一連の曲線2910を図示している。各曲線2910は単一の濡れ性状態についての毛管圧力を示す走査曲線である。画像レジストレーションによって、試料のフラッディングサイクル（たとえば1次排出、1次自発吸収、1次強制吸収、2次排出、2次吸収などの後の流体飽和状態）中の流体飽和状態において生じる変化を3Dで視覚化及び定量化することが可能となる。グラフ2910は、特定の一連の実験中に起こる様々な考えられ得る経路−たとえば経路A→B→C→D−を識別している。この一連の実験でわずかな変化が起こった結果、経路がA→E→F→Gとなることもあり得る。一組の実験は以下の手順のうちの少なくとも1つを有して良い。
1. グラフ2910の状態Aでは、岩石試料は水によって十分飽和している。
2. 1に続いて、前記岩石へのオイル又はガスの注入が、特定の最終注入圧力になるまで行われる。その結果、たとえばB（高圧）又はE（低圧）のような状態となる。
3. その後、系を濡れる流体と接触させながらその系を緩和させることによって、前記オイル又はガスは「自発的に」回復する。それによってC又はFのような状態がそれぞれ得られる。
4. その後前記試料は、強制的な濡れ性流体の注入を受ける。その結果、状態D又は状態Dになることができる。
5. 当該サイクルは手順2から任意の回数繰り返されて良い。
【０１１５】
各サイクルの各段階での流体構造の差異を定量化することができる。それにより調査対象である多相流動特性を推定することが可能となる。上述したように、2910は単一の濡れ性状態について毛管圧力を示す走査曲線を示している。複数のサイクルの実験が、他の流体によって、かつ様々な条件下において、同一のコア材料上でさらに行われて良い。そのような実験が行われることで、異なる一組の毛管圧力/水飽和の走査曲線が得られる。これらの状態は撮像することが可能である。流体飽和状態の差異が撮像され、かつ定量的に比較することができる。図30-31に図示された飽和状態は、毛管圧力走査曲線2910上において示された点A-Dに基づく。図30は、状態Aと状態Bとの間での流体分布の変化を図示している。状態Aでは、試料は水3010で飽和している。状態Bは、圧力をかけた状態で、遺留水3012が飽和するまえオイル（又はガス）を試料へ注入した後のレジストレーションされた3D体積から得られた、同一の薄片である。3012では、オイルが試料へ入り込む孔が黒で、かつ残りの水の飽和は小さな（わずかに明るい）孔である。定量的なレジストレーションは、状態Bで注入された流体の量を特定する−これは体積測定に基づいて測定することができる−だけではなく、両流体の詳細な孔スケール分布を直接観察することも可能にする。2つの流体の相対流動能力は、各相における相対透水率のシミュレーションによって3D画像で定量化することができる。そのような定量的な解析は、岩石の初期の炭化水素のフローポテンシャルの推定を行うことができる。前記フローポテンシャルは重要な製造上の示唆となるパラメータであると考えられる。
【０１１６】
図31及び図32Aの状態CとDは、状態Bで観察される初期（遺留水）の飽和からのオイル又はガスの自発変位及び強制変位後にそれぞれ明確になる。ここでは、「自発的な」流体の注入と強制的な水の注入の両方によってオイルが回復した後、どの程度オイル又はガスが試料中に残っているのかを観察することが可能である。この観察は、残りのオイル及びガスの位置を識別するので非常に重要である。残りのオイル又はガスの分布に関する知見は、さらなるフラッディング/処理の候補と考えられ得るのか否かを特定する上で助けとなる。その理由は、残りのオイル及びガスは、多数の既知の候補となる手法（たとえば第三次フラッディング、低塩分フラッディング、化学フラッディング、界面活性剤の注入、ガスの注入等）によって移動可能であるためである。これらの手法による流体の移動もまた、この方法によって直接的に検討することができる。
【０１１７】
図32Aは3D断層画像から得られた同一の薄片を図示している。前記3D断層画像は、様々な濡れ性条件下フラッディングを行った結果得られた様々な残留炭化水素分布を表している。3210は、水が濡れる条件下で水：オイル系へ水をフラッディングさせた後の残留炭化水素分布を示している。3212は、混合物の濡れ性(MW)条件下で水：オイル系へ水をフラッディングさせた後の残留炭化水素分布を示している。流体の分布はそれぞれ異なり、かつ定量化可能である。たとえば図32Bのグラフ3214乃至3220は、これらの様々な濡れ性条件下での残留炭化水素の小塊のサイズ（水のフラッディング後でのオイルを標的とする上で重要である）変化を定量化することができる。図32Bでは、残留炭化水素の小塊のサイズは、自発的吸水(SI)及び強制吸水(FI)の後に、水が濡れる条件（グラフ3218）から混合物が濡れる条件（グラフ3220）で変化する。
【０１１８】
図33は、どのようにして3D画像レジストレーションが、多数の長さスケールでの試料の特性を特定するのに利用可能であるのかを図示している。多くの材料は、（nmからcmまで）様々な寸法の3D構造部位を含んでいる。一の単一3D画像では、複数スケールの全てで全情報を捉えることができない（たとえばそのような画像は100000003すなわち200兆バイト（200テラバイト）を必要とする）。多数の長さスケールを探索するのに必要なギャップの架橋を助けるため、撮像は複数のスケールで行われて良く、かつ、定量的なレジストレーションは、大きな領域上での微細構造のマッピングに用いられて良い。
【０１１９】
mmからnmスケールにまでわたる部位を有する試料の画像の例が図33に図示されている。薄片3310によって表される3D画像が、20μmのボクセルサイズにて4cm直径の視野で得られる。これにより、20μmの材料部位と4cmの材料部位との識別が可能となる。しかし高解像度で部位を探索するため、4cmのスケールで本来撮像された試料の1つ以上の部分集合が再撮像されなければならない。この画像は本来の画像3310の部分集合で、かつ2.5μmの解像度で得られた。ここで部位は、改善された解像度で識別可能である。前記改善された解像度は元の解像度の約8倍よりも良好である。これら2つの画像から得られた情報は、3Dの画像データの組を得ることによって合成されて良い。前記3Dの画像データの組は、画像3310の各ボクセルに対してマッピングされた画像3312の高解像度画像データを有する。
【０１２０】
図34は、大きな画像3410へ埋め込まれる高解像度画像3412を図示している。高解像度画像3412は、画像3410内での部位の識別を容易にする。さもなければこれは、多数の複製後での「解像」が難しくなると思われる。
【０１２１】
図35は、各異なる解像度を有する複数の画像の同一の再分割した領域を並べて比較したものを図示している。低解像度の画像3510では解像されない部位が高解像度画像3512では識別することができる。画像レジストレーションの厳密かつ定量的な性質を有するため、高解像度の画像で識別された部位を利用することで、大きな画像中の不十分すなわち解像できない部位を「加える」ことができる。画像2510中で識別されるグレースケールは、高解像度画像3512内で明確に識別される相へ直接マッピングされて良い。たとえば画像3512内で明確に観察される有孔率と孔のサイズは、左側の画像中で識別されるグレースケール値に相関付けられて良い。図23及び図24に示された3D体積への2Dデータの伝播/合成と同様に、3412と3512から得られる関連する3D構造の情報（たとえば孔のサイズ、孔の形状、グレイン構造）は、より大きな体積（それそれ3410と3510）で画像データ全体にわたって伝播されて良い。
【０１２２】
F. 他の態様
開示された方法及びシステムは、3DのX線マイクロCTデータへの、（鉱物学及び探索のサブミクロンスケールである）従来の顕微鏡手法から取得可能な詳細な構造の情報の結合を助ける。特に、開示された方法は、マイクロCTから得られた3D画像への詳細な2D顕微鏡画像の最適なレジストレーションを可能にする。この結果、従来のマイクロCTイメージングの能力は顕著に拡張する。
【０１２３】
それに加えて、様々な濡れ性条件及び飽和状態下における同一の3D試料上での3D孔スケール分布を測定する方法並びにシステムが開示されている。当該方法は、ビームラインから有孔性試料を除去する手順、実験室内で系についての一組の実験を行う手順、及び、その後前記試料を再撮像することで、3D画像のボクセル−ボクセルのレジストレーションを行う手順を可能にする。
【０１２４】
画像データを処理する開示された方法を適用する結果、試料の1つ以上の特性に関する情報を有する合成画像データが得られる。得られた情報は、各独立した画像から得ることのできる情報に追加される。さらに情報は画像間の重なり領域を超えて補間されて良いので、当該方法は、レジストレーション処理に関与する本来の画像から得ることのできる結合された情報に追加される情報を得ることを可能にする。
【０１２５】
さらに当該方法は、試料の撮像領域全体にわたる−ひいては試料全体にわたる−その試料の特性の空間分布を推定することを可能にする。ここで「特性」という語は、非常に広い意味で用いられていて、かつ試料の情報を有する如何なる特徴をも含むと解される。そのような特性の例は以下の1つ以上を含んで良い。それは、物理的特性、構造特性、化学特性、鉱物学的特性、透水率、音響特性、グレインサイズ、グレイン境界、セメント相、地質学上外見の種類と地質学上の岩石の種類、鉱物相の構造、又は試料の孔への流体分布である。特に開示された方法のこのような有利な能力は、2つの3D画像の処理のレジストレーションの場合にも適用可能である。
【特許請求の範囲】
【請求項１】
試料の画像データを処理する方法であって：
前記試料の少なくとも一部が重なっている複数の空間領域の第1画像と第2画像とのレジストレーションを行う手順；及び、
前記のレジストレーションされた画像からのデータを処理することで、前記第1画像と第2画像から得られる画像に加えられる前記試料についての情報を有する合成画像データを取得する手順；
を有する方法。
【請求項２】
前記情報が、前記試料の特性の空間分布の定量的な推定を有する、請求項1に記載の方法。
【請求項３】
前記レジストレーションが、漸次解像度が高くなる画像に基づいて変換パラメータを計算する繰り返し処理を実行する手順を有する、請求項2に記載の方法。
【請求項４】
前記第1画像が3D画像で、かつ前記第2画像が相対的に高い解像度の2D画像で、
さらに、前記レジストレーションは、前記2D画像データの、前記3D画像データの対応する薄片に対する座標の位置合わせを有する、
請求項1乃至3のうちいずれか一項に記載の方法。
【請求項５】
さらに、前記薄片についての合成画像データが、前記2D画像に係るピクセルデータをマッピングすることで、前記3D画像の薄片の対応するレジストレーションされたボクセルを追加することによって得られる、請求項4に記載の方法。
【請求項６】
さらに、前記2D画像と前記薄片から得られた画像データが、前記3D画像によって撮像された少なくとも空間領域全体にわたる前記特性の空間分布についての定量的な推定が得られるように処理される、請求項5に記載の方法。
【請求項７】
前記2D画像データが、ある範囲の長さスケールにわたって前記定量的な推定を供するように前処理される、請求項4乃至6のうちいずれか一項に記載の方法。
【請求項８】
前記2D画像が、前記3D画像に対してレジストレーションされる前に湾曲歪み補正される、請求項4乃至7のうちいずれか一項に記載の方法。
【請求項９】
さらに、前記第1画像及び前記第2画像が3D画像である、請求項1乃至3のうちいずれか一項に記載の方法。
【請求項１０】
さらに、前記第1画像が前記第2画像よりも高い解像度で、かつ前記合成画像データが、前記第1の3D画像に係るピクセル画像データをマッピングすることで、前記第2の3D画像の対応するレジストレーションされたボクセルを追加することによって得られる、請求項9に記載の方法。
【請求項１１】
さらに、前記第1の3D画像が、前記第2の3D画像によって撮像される少なくとも空間領域全体にわたる前記特性の空間分布についての定量的な推定が得られるように処理される、請求項10に記載の方法。
【請求項１２】
さらに、前記第1の3D画像と前記第2の3D画像が異なる時点で撮像され、かつ
前記の取得された合成画像データが第3の3D画像を生成し、
前記第3画像は、前記第1画像の画像データと前記第2画像の画像データとの差異を示すデータによって追加される、
請求項9に記載の方法。
【請求項１３】
前記特性が、物理的特性、構造特性、化学特性、鉱物学的特性、透水率、音響特性、グレインサイズ、グレイン境界、セメント相、地質学上外見の種類と地質学上の岩石の種類、鉱物相の構造、又は試料の孔への流体分布のうちの少なくとも1つを有する、請求項1乃至12のうちいずれか一項に記載の方法。
【請求項１４】
前記合成画像データを用いることで前記の撮像された試料の物理的特性のシミュレーションを行う手順をさらに有する請求項1乃至13のうちいずれか一項に記載の方法であって、
前記物理的特性は、前記試料材料の様々な流体相の相対透水率、誘電応答、NMR応答、及び弾性のうちの少なくとも1つを有する、
方法。
【請求項１５】
前記3D画像がマイクロCT撮像によって取得され、かつ
前記2D画像が以下の手法のうちのいずれか一によって取得され、
前記以下の手法とは：
走査型電子顕微鏡(SEM)；
光学顕微鏡；
集束イオンビームSEM；
レーザー共焦点顕微鏡；
エネルギー分散X線解析；
X線蛍光；
X線光電子分光；
2次イオン質量分析分光；
赤外及びラマン分光；
UV及び可視分光；
走査型音響顕微鏡；並びに、
原子間力顕微鏡；
である、
請求項4乃至8のうちいずれか一項に記載の方法。
【請求項１６】
前記画像レジストレーションが、多解像度多開始点の大域最適化法を用いることによって実行される、請求項1乃至15のうちいずれか一項に記載の方法。
【請求項１７】
前記最適化法を実装する計算計画が、タスクの並列化及びデータの並列化のうちの少なくとも1つを有する、請求項16に記載の方法。
【請求項１８】
前記画像レジストレーションの計算ルーチンが、分配されたメモリ並列アルゴリズムを実装する、請求項17に記載の方法。
【請求項１９】
前記第1画像及び前記第2画像が、相互に対して任意に配向し、並進し、かつ/又は伸縮する、請求項1乃至18のうちいずれか一項に記載の方法。
【請求項２０】
さらに、前記の取得された合成画像データが前記試料に関する情報を有し、
該情報は、前記第1画像及び前記第2画像が独立するものとみなされるときには、前記第1画像及び前記第2画像から取得可能な情報に追加される、
請求項1乃至19のうちいずれか一項に記載の方法。
【請求項２１】
請求項1乃至20のうちいずれか一項に記載の方法を実行することが可能なコンピュータプログラム。
【請求項２２】
請求項20に記載のコンピュータプログラムを有するコンピュータによる読み取りが可能な媒体を有するコンピュータプログラム製品。
【請求項２３】
試料の画像データを処理する電子システムであって、
当該システムは：
前記試料の少なくとも部分的に重なる複数の空間領域の第1画像と第2画像を取得するように備えられた少なくとも1つの撮像装置；
前記の取得された画像を記憶するように備えられた記憶デバイス；及び、
前記の取得された画像のデータを取得して、レジストレーションすることで合成画像データを得るように備えられた少なくとも1つのマイクロプロセッサ；
を有し、
前記合成画像データは前記試料に関する情報を有し、
前記情報は前記第1画像及び前記第2画像から取得可能な情報に追加される、
システム。
【請求項２４】
前記少なくとも1つの撮像装置が3D撮像装置を有する、請求項23に記載のシステム。
【請求項２５】
前記少なくとも1つの撮像装置が2D撮像装置及び3D撮像装置を有し、
前記2D撮像装置は、前記3D撮像装置によって得られる解像度よりも実質的に高い解像度の画像を生成することが可能である、
請求項23に記載のシステム。
【請求項１】
試料の画像データを処理する方法であって：
前記試料の少なくとも一部が重なっている複数の空間領域の第1画像と第2画像とのレジストレーションを行う手順；及び、
前記のレジストレーションされた画像からのデータを処理することで、前記第1画像と第2画像から得られる画像に加えられる前記試料についての情報を有する合成画像データを取得する手順；
を有する方法。
【請求項２】
前記情報が、前記試料の特性の空間分布の定量的な推定を有する、請求項1に記載の方法。
【請求項３】
前記レジストレーションが、漸次解像度が高くなる画像に基づいて変換パラメータを計算する繰り返し処理を実行する手順を有する、請求項2に記載の方法。
【請求項４】
前記第1画像が3D画像で、かつ前記第2画像が相対的に高い解像度の2D画像で、
さらに、前記レジストレーションは、前記2D画像データの、前記3D画像データの対応する薄片に対する座標の位置合わせを有する、
請求項1乃至3のうちいずれか一項に記載の方法。
【請求項５】
さらに、前記薄片についての合成画像データが、前記2D画像に係るピクセルデータをマッピングすることで、前記3D画像の薄片の対応するレジストレーションされたボクセルを追加することによって得られる、請求項4に記載の方法。
【請求項６】
さらに、前記2D画像と前記薄片から得られた画像データが、前記3D画像によって撮像された少なくとも空間領域全体にわたる前記特性の空間分布についての定量的な推定が得られるように処理される、請求項5に記載の方法。
【請求項７】
前記2D画像データが、ある範囲の長さスケールにわたって前記定量的な推定を供するように前処理される、請求項4乃至6のうちいずれか一項に記載の方法。
【請求項８】
前記2D画像が、前記3D画像に対してレジストレーションされる前に湾曲歪み補正される、請求項4乃至7のうちいずれか一項に記載の方法。
【請求項９】
さらに、前記第1画像及び前記第2画像が3D画像である、請求項1乃至3のうちいずれか一項に記載の方法。
【請求項１０】
さらに、前記第1画像が前記第2画像よりも高い解像度で、かつ前記合成画像データが、前記第1の3D画像に係るピクセル画像データをマッピングすることで、前記第2の3D画像の対応するレジストレーションされたボクセルを追加することによって得られる、請求項9に記載の方法。
【請求項１１】
さらに、前記第1の3D画像が、前記第2の3D画像によって撮像される少なくとも空間領域全体にわたる前記特性の空間分布についての定量的な推定が得られるように処理される、請求項10に記載の方法。
【請求項１２】
さらに、前記第1の3D画像と前記第2の3D画像が異なる時点で撮像され、かつ
前記の取得された合成画像データが第3の3D画像を生成し、
前記第3画像は、前記第1画像の画像データと前記第2画像の画像データとの差異を示すデータによって追加される、
請求項9に記載の方法。
【請求項１３】
前記特性が、物理的特性、構造特性、化学特性、鉱物学的特性、透水率、音響特性、グレインサイズ、グレイン境界、セメント相、地質学上外見の種類と地質学上の岩石の種類、鉱物相の構造、又は試料の孔への流体分布のうちの少なくとも1つを有する、請求項1乃至12のうちいずれか一項に記載の方法。
【請求項１４】
前記合成画像データを用いることで前記の撮像された試料の物理的特性のシミュレーションを行う手順をさらに有する請求項1乃至13のうちいずれか一項に記載の方法であって、
前記物理的特性は、前記試料材料の様々な流体相の相対透水率、誘電応答、NMR応答、及び弾性のうちの少なくとも1つを有する、
方法。
【請求項１５】
前記3D画像がマイクロCT撮像によって取得され、かつ
前記2D画像が以下の手法のうちのいずれか一によって取得され、
前記以下の手法とは：
走査型電子顕微鏡(SEM)；
光学顕微鏡；
集束イオンビームSEM；
レーザー共焦点顕微鏡；
エネルギー分散X線解析；
X線蛍光；
X線光電子分光；
2次イオン質量分析分光；
赤外及びラマン分光；
UV及び可視分光；
走査型音響顕微鏡；並びに、
原子間力顕微鏡；
である、
請求項4乃至8のうちいずれか一項に記載の方法。
【請求項１６】
前記画像レジストレーションが、多解像度多開始点の大域最適化法を用いることによって実行される、請求項1乃至15のうちいずれか一項に記載の方法。
【請求項１７】
前記最適化法を実装する計算計画が、タスクの並列化及びデータの並列化のうちの少なくとも1つを有する、請求項16に記載の方法。
【請求項１８】
前記画像レジストレーションの計算ルーチンが、分配されたメモリ並列アルゴリズムを実装する、請求項17に記載の方法。
【請求項１９】
前記第1画像及び前記第2画像が、相互に対して任意に配向し、並進し、かつ/又は伸縮する、請求項1乃至18のうちいずれか一項に記載の方法。
【請求項２０】
さらに、前記の取得された合成画像データが前記試料に関する情報を有し、
該情報は、前記第1画像及び前記第2画像が独立するものとみなされるときには、前記第1画像及び前記第2画像から取得可能な情報に追加される、
請求項1乃至19のうちいずれか一項に記載の方法。
【請求項２１】
請求項1乃至20のうちいずれか一項に記載の方法を実行することが可能なコンピュータプログラム。
【請求項２２】
請求項20に記載のコンピュータプログラムを有するコンピュータによる読み取りが可能な媒体を有するコンピュータプログラム製品。
【請求項２３】
試料の画像データを処理する電子システムであって、
当該システムは：
前記試料の少なくとも部分的に重なる複数の空間領域の第1画像と第2画像を取得するように備えられた少なくとも1つの撮像装置；
前記の取得された画像を記憶するように備えられた記憶デバイス；及び、
前記の取得された画像のデータを取得して、レジストレーションすることで合成画像データを得るように備えられた少なくとも1つのマイクロプロセッサ；
を有し、
前記合成画像データは前記試料に関する情報を有し、
前記情報は前記第1画像及び前記第2画像から取得可能な情報に追加される、
システム。
【請求項２４】
前記少なくとも1つの撮像装置が3D撮像装置を有する、請求項23に記載のシステム。
【請求項２５】
前記少なくとも1つの撮像装置が2D撮像装置及び3D撮像装置を有し、
前記2D撮像装置は、前記3D撮像装置によって得られる解像度よりも実質的に高い解像度の画像を生成することが可能である、
請求項23に記載のシステム。
【図１】
【図２Ａ】
【図２Ｂ】
【図３】
【図４Ａ】
【図４Ｂ】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１Ａ】
【図１１Ｂ】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６Ａ】
【図１６Ｂ】
【図１７】
【図１８】
【図１９】
【図２０】
【図２１】
【図２２】
【図２３】
【図２４Ａ】
【図２４Ｂ】
【図２５】
【図２６】
【図２７】
【図２８】
【図２９】
【図３０】
【図３１】
【図３２Ａ】
【図３２Ｂ】
【図３３】
【図３４】
【図３５】
【図３６Ａ】
【図３６Ｂ】
【図２Ａ】
【図２Ｂ】
【図３】
【図４Ａ】
【図４Ｂ】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１Ａ】
【図１１Ｂ】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６Ａ】
【図１６Ｂ】
【図１７】
【図１８】
【図１９】
【図２０】
【図２１】
【図２２】
【図２３】
【図２４Ａ】
【図２４Ｂ】
【図２５】
【図２６】
【図２７】
【図２８】
【図２９】
【図３０】
【図３１】
【図３２Ａ】
【図３２Ｂ】
【図３３】
【図３４】
【図３５】
【図３６Ａ】
【図３６Ｂ】
【公表番号】特表２０１１−５２５２７１（Ｐ２０１１−５２５２７１Ａ）
【公表日】平成２３年９月１５日（２０１１．９．１５）
【国際特許分類】
【出願番号】特願２０１１−５０９８２０（Ｐ２０１１−５０９８２０）
【出願日】平成２１年５月２２日（２００９．５．２２）
【国際出願番号】ＰＣＴ／ＡＵ２００９／０００６４１
【国際公開番号】ＷＯ２００９／１４０７３８
【国際公開日】平成２１年１１月２６日（２００９．１１．２６）
【公序良俗違反の表示】
（特許庁注：以下のものは登録商標）
１．フロッピー
【出願人】（５０６２７０７５４）ザ オーストラリアン ナショナル ユニヴァーシティ (2)
【Ｆターム（参考）】
【公表日】平成２３年９月１５日（２０１１．９．１５）
【国際特許分類】
【出願日】平成２１年５月２２日（２００９．５．２２）
【国際出願番号】ＰＣＴ／ＡＵ２００９／０００６４１
【国際公開番号】ＷＯ２００９／１４０７３８
【国際公開日】平成２１年１１月２６日（２００９．１１．２６）
【公序良俗違反の表示】
（特許庁注：以下のものは登録商標）
１．フロッピー
【出願人】（５０６２７０７５４）ザ オーストラリアン ナショナル ユニヴァーシティ (2)
【Ｆターム（参考）】
[ Back to top ]