https://www.academia.edu/download/74363184/bf0241105020211108-2478-5qo6ho.pdf 要旨 三次元光弾性体は、ある角度 θ における線形リターダー δ と純回転子 χ から構成される光学的に等価なモデルとして表すことができる。これらは、それぞれ特性遅相差 δ、一次および二次の特性方向 θ および θ+χ と呼ばれる。これまで、これらの特性パラメータは、手動による逐点的なデータ収集法でしか求められず、非常に手間と時間を要していた。 そこで、三次元あるいは統合光弾性における3つの特性パラメータを決定するための、自動化された位相シフト法が開発された。 6枚の位相シフト画像から δ、θ、そして θ+χ を求めるための数式が導出されている。これらの画像はCCDカメラを用いて収集され、標準的なパーソナルコンピュータで全視野データを処理する。この新しい方法により、比較的短時間で3つの特性パラメータすべての正確な全視野マップを得ることが可能になり、光弾性データの全視野トモグラフィー再構成が現実的に実現できるようになった。 キーワード ― 特性パラメータ、統合光弾性、位相シフト法、ポアンカレ等価定理  

序論（Introduction 日本語訳）

一般的な三次元光弾性体では、主応力差の大きさと方向は光路に沿って連続的に変化する。この光弾性体は、ある角度 θ における遅相 δ を持つ線形リターダーと、純粋な回転子 χ から構成される光学的に等価なモデルで表すことができる。Aben² はこれらを三次元光弾性体の「特性パラメータ」と呼び、δ を特性遅相差、θ および (θ+χ) をそれぞれ一次および二次の特性方向とした。 これまで、これらの特性パラメータは手動による逐点的な測定法によってしか決定できず、煩雑かつ時間を要していた。たとえば、2つの特性方向を決定する過程では、透過型ポラリスコープ内の光学素子を回転させて、注目する各点で光強度の最小値と最大値を得る必要がある。特性遅相差は二次元光弾性で用いられる Tardy 法で決定可能であるが、それでも逐点測定が必要であり、解析は限定的で時間もかかる。 本研究の目的は、二次元光弾性で使用されてきた自動化位相シフト法を拡張し、三次元あるいは統合光弾性における3つの特性パラメータを決定できるようにすることで、これらの問題を克服することである。これは Patterson と Wang⁶ によって開発された技術を発展させるものである。彼らは CCD カメラで位相シフト画像を取得し、それを標準的なパーソナルコンピュータで処理することで、等色線次数や等傾線角度の全視野マップを得ていた。本研究ではこの手法をさらに改良し、比較的短時間で3つの特性パラメータすべての正確な全視野マップを取得できるようにした。これにより光弾性トモグラフィの実現が現実味を帯びてきた。 三次元光弾性モデルの複数の方向から δ、θ、χ を取得すれば、応力凍結モデルを切断することなくトモグラフィ再構成によって主応力差とその方向を求めることが可能となる。したがって、モデルは破壊されず、切断過程でデータを失うこともない。このように位相シフト法を用いるアイデアは、著者らによって 1998 年に初めて提案され⁷⁸、翌 1999 年には類似手法が報告されている⁹。 さらに、この手法は光弾性トモグラフィへの応用に限らず、特性パラメータの決定を必要とする散乱光光弾性¹⁰、磁気光弾性¹¹、複屈折光ファイバーセンシング¹² ¹³ などの手法を強化し、三次元複屈折流れの測定¹⁴ を可能にする。 本論文では、位相シフト式の理論展開を概説し、新しい手法で処理した3つのデータ例を提示する。  

特性パラメータの決定

位相シフト式

図1に示すように、任意の方向に配置された1/4波長板と出力偏光子（アナライザ）を備えたポラリスコープを考える。光源からの非偏光光の入力はストークスベクトル S であり、 S=(1000)(1)S = \begin{pmatrix} I \\ 0 \\ 0 \\ 0 \end{pmatrix} \tag{1} で表される。 垂直基準軸に対して角度 φ をなす理想的線形偏光子の Mueller 行列 は次式で与えられる。 Pϕ=12(1cos⁡2ϕsin⁡2ϕ0cos⁡2ϕcos⁡22ϕcos⁡2ϕsin⁡2ϕ0sin⁡2ϕcos⁡2ϕsin⁡2ϕsin⁡22ϕ00000)(2)P_\phi = \frac{1}{2} \begin{pmatrix} 1 & \cos 2\phi & \sin 2\phi & 0 \\ \cos 2\phi & \cos^2 2\phi & \cos 2\phi \sin 2\phi & 0 \\ \sin 2\phi & \cos 2\phi \sin 2\phi & \sin^2 2\phi & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} \tag{2} 速軸が基準軸となす角 φ の 1/4 波長板 の Mueller 行列は Qϕ=(10000cos⁡22ϕcos⁡2ϕsin⁡2ϕ−sin⁡2ϕ0cos⁡2ϕsin⁡2ϕsin⁡22ϕcos⁡2ϕ0sin⁡2ϕ−cos⁡2ϕ0)(3)Q_\phi = \begin{pmatrix} 1 & 0 & 0 & 0 \\ 0 & \cos^2 2\phi & \cos 2\phi \sin 2\phi & -\sin 2\phi \\ 0 & \cos 2\phi \sin 2\phi & \sin^2 2\phi & \cos 2\phi \\ 0 & \sin 2\phi & -\cos 2\phi & 0 \end{pmatrix} \tag{3} で表される。 三次元試料は、基準軸に対して角度 θ の速軸を持つ遅相差 δ の 線形リターダー と、回転子 χ に光学的に等価である。リターダーの Mueller 行列は Mθ(δ)=(10000cos⁡22θ+sin⁡22θcos⁡δsin⁡2θcos⁡2θ(1−cos⁡δ)−sin⁡2θsin⁡δ0sin⁡2θcos⁡2θ(1−cos⁡δ)sin⁡22θ+cos⁡22θcos⁡δcos⁡2θsin⁡δ0sin⁡2θsin⁡δ−cos⁡2θsin⁡δcos⁡δ)(4)M_\theta(\delta) = \begin{pmatrix} 1 & 0 & 0 & 0 \\ 0 & \cos^2 2\theta + \sin^2 2\theta \cos \delta & \sin 2\theta \cos 2\theta (1 - \cos \delta) & -\sin 2\theta \sin \delta \\ 0 & \sin 2\theta \cos 2\theta (1 - \cos \delta) & \sin^2 2\theta + \cos^2 2\theta \cos \delta & \cos 2\theta \sin \delta \\ 0 & \sin 2\theta \sin \delta & -\cos 2\theta \sin \delta & \cos \delta \end{pmatrix} \tag{4} 回転子の Mueller 行列は Rχ=(10000cos⁡2χ−sin⁡2χ00sin⁡2χcos⁡2χ00001)(5)R_\chi = \begin{pmatrix} 1 & 0 & 0 & 0 \\ 0 & \cos 2\chi & -\sin 2\chi & 0 \\ 0 & \sin 2\chi & \cos 2\chi & 0 \\ 0 & 0 & 0 & 1 \end{pmatrix} \tag{5} で与えられる。 したがって、ポラリスコープ出力における光のストークスベクトル S_out は次式で表される。 Sout=Pϕ3 Qϕ Rχ Mθ(δ) Q0 P0 S(6)S_{\text{out}} = P_{\phi_3} \, Q_{\phi} \, R_\chi \, M_\theta(\delta) \, Q_0 \, P_0 \, S \tag{6}

強度式

入力および出力素子の配列を表1に従って設定した場合、観測される強度は次式のようになる。 I1=a22(7a)I_1 = \frac{a^2}{2} \tag{7a} I2=a22(1+sin⁡δcos⁡2(θ+χ))(7b)I_2 = \frac{a^2}{2} \left( 1 + \sin \delta \cos 2(\theta + \chi) \right) \tag{7b} I3=a22(1−sin⁡δcos⁡2(θ+χ))(7c)I_3 = \frac{a^2}{2} \left( 1 - \sin \delta \cos 2(\theta + \chi) \right) \tag{7c} I4=a22(1−sin⁡δsin⁡2(θ+χ))(7d)I_4 = \frac{a^2}{2} \left( 1 - \sin \delta \sin 2(\theta + \chi) \right) \tag{7d} I5=a22(1+sin⁡δsin⁡2(θ+χ))(7e)I_5 = \frac{a^2}{2} \left( 1 + \sin \delta \sin 2(\theta + \chi) \right) \tag{7e} I6=a22(1+cos⁡δ)(7f)I_6 = \frac{a^2}{2} \left( 1 + \cos \delta \right) \tag{7f} ここで a は比例定数である。

特性パラメータ算出式

これら6枚の位相シフト画像を組み合わせることで、特性パラメータは次式で決定できる。 一次・二次特性方向: (θ+χ)=12tan⁡−1(I2−I3I3−I2−2I1)(8)(\theta + \chi) = \frac{1}{2} \tan^{-1} \left( \frac{I_2 - I_3}{I_3 - I_2 - 2I_1} \right) \tag{8} 特性遅相差: δ=tan⁡−1(sin⁡2(θ+χ)(I2−I3)+I3−2I4sin⁡2θcos⁡2(θ+χ))(9)\delta = \tan^{-1} \left( \frac{\sin^2(\theta + \chi) (I_2 - I_3) + I_3 - 2I_4}{\sin 2\theta \cos 2(\theta + \chi)} \right) \tag{9} 一次特性方向 θ: θ=12tan⁡−1(sin⁡2(θ+χ)sin⁡δ1−2I6I2+I3−cos⁡2(θ+χ)sin⁡2θcos⁡δ)(10)\theta = \frac{1}{2} \tan^{-1} \left( \frac{ \sin 2(\theta + \chi) \sin \delta }{ 1 - \frac{2I_6}{I_2 + I_3} - \cos 2(\theta + \chi) \sin 2\theta \cos \delta } \right) \tag{10} 補助式: sin⁡2θ=(2I5I2+I3)−1,sin⁡δ=2I5I2+I3−1(11)\sin 2\theta = \left( \frac{2I_5}{I_2 + I_3} \right) - 1, \quad \sin \delta = \frac{2I_5}{I_2 + I_3} - 1 \tag{11}

表1 — 入出力1/4波長板とアナライザの配置

画像番号 Qp Qφ Pφ3
1 π/4 π/4 π/4
2 π/4 0 0
3 π/4 0 π/2
4 π/4 0 π/4
5 π/2 0 π/4
6 π/2 0 π/2
  これらの式から、理論的には 5枚の位相シフト画像（I₁～I₅） だけで3つの特性パラメータを得ることが可能であることがわかる。 しかし、式 (8)〜(10) の解は周期的なパラメータのマップであり、アンラッピングアルゴリズムを用いて連続的なマップに変換する必要がある。 Wang と Patterson¹⁶ によって記述されたアルゴリズムでは arctan（逆正接関数） が使用されており、これによってデータマップ中の転換点を容易に特定できる。 一方、もし5枚の画像のみを使用する場合には arccos（逆余弦関数） が必要となり、転換点の識別は容易ではない。 したがって、本研究では6枚の画像を用いてこれらのアルゴリズムを実装することにした。 類似の方法に関する最近の論文において、Mangal と Ramesh⁹ は「統合光弾性に位相シフト法を適用するには12枚の画像が必要である」と述べているが、これは明らかに当てはまらない。  

復調とアンラップ（Demodulation and Unwrapping）

特性パラメータのアンラップ処理では、互いに干渉し合うため複数の問題が生じ、またデータ中に未定義領域が現れる。これらは、圧縮下の円板から得られたデータを示す図2で説明できる。最初の2つの問題は、Patterson と Wang によって二次元光弾性用に開発された位相シフト手法でも見られたものである。 (θ+χ) の周期性は、相対遅相差とその一次導関数の符号を反転させ、特に二次特性方向 (θ+χ) = ±45° の位置で顕著となる（図2(a)参照）。また、図2(a)から、相対遅相差 δ がゼロのときには二次特性方向が未定義になることも確認できる。これらの領域は容易に特定でき、Mawatari ら¹⁷ が提案したようにデータを外挿して補間することができる。 図2(b)は3つの特性パラメータに関する周期データを示している。未定義領域が一次特性方向にも存在することが確認できる。さらに δ の周期性は、θ の符号およびその一次導関数の反転を引き起こす。これは δ の周期端で発生し、この不連続性は遅相差の情報から特定できる。ただし、この処理は (θ+χ) によって生じる δ の段差を取り除いた後に行う必要がある。 位相データの復調は Wang と Patterson¹⁶ による手法を用いて実施される。θ の正しい復調のためには、二次特性方向 (θ+χ) と同相である必要がある。θ の復調とアンラップの結果は図2(c) および図2(d)に示されている。図2(c)には周期的な θ と復調後の θ が、図2(d)にはアンラップされた θ データとアンラップされた (θ+χ) の比較が示されている。 したがって、位相データの復調とアンラップの手順は以下の通りである：
  1. 周期的な (θ+χ)、δ、θ を求める
  2. (θ+χ) を用いて δ を復調する
  3. (θ+χ) と δ をアンラップする
  4. 復調された周期的 δ を用いて θ を復調する
  5. θ をアンラップする

実験手順（Experimental Procedure）

上述の理論は以下の実験手順に実装された。ポラリスコープの配置は図3に示されており、データは CCD カメラを用いて収集され、フレームグラバー付きのパソコンに接続された。三次元複屈折モデルは、モデルと屈折率を一致させた浸液を満たした透明タンク内の標準透過型ポラリスコープに設置される。モデルは応力凍結状態または実荷重状態とすることができる。入力光は初期状態で円偏光化されている。 6枚の画像は順に、表1に示された光学素子の配置に回転させながら取得され、各画像はフレームグラバーを用いて記録される。この結果、256×256ピクセルのデジタルデータマップが得られる。これらの画像は式(8)〜(10)を用いて処理され、アンラップされることで3つの特性パラメータの全視野マップが得られる。

例（Examples）

圧縮円板（Disc in Compression）

アルゴリズムの最初のテストは、数学ソフトウェア（MATHCAD Version Plus 7.0 Professional Edition）を用いたシミュレーション実験で行われた。直径方向に圧縮された光弾性円板が生成され、表1に示されたポラリスコープ素子の配置に基づき6枚の位相シフト画像が作成された。これらの画像は、式(8)〜(10)を解き、3つの特性パラメータの全視野マップを生成する専用プログラムに入力された。 図4(a)には円板の1/4領域における周期的マップが、図4(b)には円板全体のアンラップされた連続マップが示されている。アンラップ後の一次特性方向と二次特性方向が一致していることが確認でき、これは二次元モデルにおいて期待される結果である。図4(b)の結果と、このモデルに対する理論解との差は計算され、図4(c)に示されている。

三次元シミュレーション（Three-dimensional Simulation）

次に三次元複屈折モデルをシミュレーションした。これは、直径方向に圧縮された3枚の光弾性円板を円筒状に組み合わせたものである（図5参照）。各円板は、荷重点が水平基準軸に対して0°、22.5°、45°となるように配置され、図5に示す方向から見たときに三次元応力場を生じるようになっている。前と同じ手順で位相シフト画像を生成し、周期的な全視野マップを解いた結果が図6に示されている。

引張板中のコーナークラック（Corner Crack in Tensile Plate）

エポキシ樹脂製の板（寸法：220×48×6.5 mm）に、直径32 mm・厚さ0.2 mmのダイヤモンドチップ丸鋸を用いてコーナークラックを加工した。鋸刃は板の表面からオフセットして切り込みを入れることで、板厚方向に沿ってクラック長が変化するようにした（図7参照）。 モデルは引張荷重を受け、応力凍結サイクルを施された。この結果、板厚を通じてクラック周囲に不均一な応力場が発生した。 板は透過型ポラリスコープ（図3）に設置され、表1に示された光学素子の6つの配置に対してナトリウム光を用いてCCDカメラで画像が記録された。位相シフト法により3つの特性パラメータの全視野マップが生成され、その結果が図8に示されている。 クラック先端付近やクラック面に隣接する領域のデータにはマスクが施された。その理由は複数ある。第一に、アンラップ処理は物理的不連続（この場合は加工されたクラック）を跨いで行えないため、この領域を除去する必要がある。第二に、この領域のデータは過剰なノイズを生じ、アンラップ処理に重大な誤差を引き起こすことがわかった。第三に、正確なアンラップには1フリンジあたり10ピクセル以上の解像度が必要であり、クラック周辺の急峻なフリンジ勾配は除去しなければならない。 以前のクラック先端研究で用いられた標準マスク（例：文献18）が適用された。これはクラック先端を中心とした半径 10ρ の円を基に構築されており、ここで ρ はクラックまたは切欠きの先端半径である。この円の内部は塑性域を含むためマスクで除去される。残りのマスク境界は、クラック経路に対して±45°の接線で定義される。得られた連続的な特性遅相差マップは、図8(b) に示すように、クラック周囲の等色線フリンジ写真と比較することができる。

考察（Discussion）

この新しい方法の主な利点は、比較的短時間で全視野データを取得できる点にある。図4に示された解像度のマップは、従来の逐点法を用いれば取得に数日を要したであろうが、現在では3つの特性パラメータすべてのマップを10分以内で得ることができる。データ収集速度の制約は、位相シフトごとに透過型ポラリスコープの光学素子を回転させる必要があることのみであり、CCDカメラは毎秒25フレームという比較的高速でデータを記録するため、この制約はほとんど無視できる。素子の回転はステッピングモータで自動化することも可能である。 全視野マップのデータ処理自体は、Pentium級のパソコンを用いればおよそ1分で完了する。性能の低いコンピュータでも実行可能であるが、その場合は処理時間が長くなる。装置は非常にシンプルであり、標準的な透過型ポラリスコープ、CCDカメラ、フレームグラバー、パソコンのみで構成される。データマップの解像度はCCDカメラのレンズ仕様にのみ依存するため、必要に応じて非常に詳細なマップを得ることが可能である。 二次元モデルを用いたアルゴリズムの検証実験は成功であり、図4の位相マップに示されている。一次特性パラメータの周期マップを全視野で復調した効果は、図4(a)(iii) と (iv) を比較することで明確である。図4(a)(iii) では θ の周期マップが示されており、δ の周期端で不連続性がはっきりと見て取れる。一方、図4(a)(iv) では不連続性の位置が特定でき、全視野マップが復調されている。周期マップはアンラップ処理により連続関数へと変換され、図4(b) に示されている。モデルが二次元であるため、回転項 χ はゼロとなり、(θ+χ) と θ のマップは一致する。図4(c) では手法の精度が示されており、特性遅相差は接触点付近を除き 0.1 フリンジ以内の精度で求められている。接触点ではフリンジ勾配が急峻で、アンラップ解像度（10ピクセル/フリンジ）を超えるため、関心領域を縮小し倍率を上げることが唯一の解決策となる。 特性方向は、未定義領域（データ中のノイズが見られる部分）や δ の周期端（図中の白色部分）を除けば、5°以内の精度で求められている。未定義領域の補間も試みられたが、接触点付近のように領域が密集している場合には限界がある。補間はノイズを取り除く効果とデータマップの真の形状を保持する効果のトレードオフを伴うため、未定義データをすべて除去することはできなかった。未定義データを完全に削除すると、元のデータが不足し、補間過程でマップの形状が歪む。別の方法として、Mangal と Ramesh が行ったようにノイズの多いデータを完全に除去する方法もあるが、この場合マップに空白や不連続が生じ、一次特性方向のアンラップ処理に誤差を生じることもある。しかし、一次・二次特性方向は同相であるべきため、両者を直接比較することで誤差の有無は明確になる。周期データのアンラップには複数の方法があるが、ここでは割愛する。 三次元応力場シミュレーションの周期的特性パラメータマップは図6に示されている。マップの解像度は非常に詳細で、未定義領域も確認できる。図6(iv) に示されるように、一次特性パラメータは適切に復調されている。 位相シフト法をコーナークラックを有するモデルに適用した最終結果は図8に示されている。図8(b)(i) の特性遅相差の全視野マップは、図8(b)(ii) のクラック先端周囲の等色線フリンジ写真とよく一致している。アンラップ処理の際、マップ上部にわずかな不連続性が生じているが、全体として本研究で提示した例は、三次元応力場の3つの特性パラメータの全視野マップを得ることが可能であることを示している。 本手法の速度と解像度は、磁気光弾性¹¹ や偏光光学時間領域反射法（POTDR）¹³ など、特性パラメータの決定を必要とするさまざまな手法に大きな可能性を与える。しかし、この手法の真の可能性は光弾性トモグラフィにあると考えられる。Aben²⁰ が提案したような全視野版の手法を用いれば、応力凍結モデルを浸液槽に設置し、光源と検出器をモデル周囲に回転させることができる。本研究の新しい位相シフト法を使えば、モデル自体を光路内で回転させることもできる。モデルは応力凍結状態でも実荷重状態でもよいが、後者の場合は荷重装置が光路を遮らないようにする必要がある。 Keshaven⁵ は最近、三次元構造部材の応力場を決定するために、複屈折モデルに漸増的な複数荷重を加え、得られる特性パラメータを利用する代替的手法を提案している。この方法の全視野版を本研究の位相シフト手法で実現するための実験も進行中である。ただし、詳細な応力パターンを得るには多数の荷重条件が必要であり、制御が難しい可能性がある。 二次元および三次元の流動複屈折については、透過型ポラリスコープや散乱光ポラリスコープを改良した装置を用い、従来は逐点的な定量解析によって研究されてきた。本研究で提案する位相シフト法は、流動複屈折に対しても全視野定量解析の可能性を与え、さらにトモグラフィ手法を実装できる可能性を持っている。

結論（Conclusions）

位相シフト法を用いて全視野の特性パラメータを取得する新しい方法を開発し、この手順で処理したデータ例を示した。本研究により、特性パラメータの周期的マップを正確に生成できること、そして連続関数を得るためのアンラップ処理においてのみ大きな誤差が生じうることが明らかとなった。

この手法によって、応力凍結モデルを切断することなく三次元応力場を決定することが現実的に可能となった。これは、多重荷重法やトモグラフィ的視点の再構成によって実現できる。

さらに、この手法は光弾性トモグラフィへの応用にとどまらず、散乱光光弾性、磁気光弾性、複屈折光ファイバーセンシングの能力を高め、三次元流動複屈折の測定を可能にする潜在力を有している。

謝辞（Acknowledgments）

著者らは、本研究に対する 英国工学・物理科学研究評議会（EPSRC） の財政的支援、および Z. E. Wang 博士 の協力に心から感謝の意を表する。